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1 Introduction 



This paper belongs to a series devoted to numerical evaluation of the multi-loop, multi-leg Feynman 
diagrams that appear in any renormalizablc quantum field theory. In pp (hereafter I) the general strategy 
has been designed and in [2] (hereafter II) a complete list of results has been derived for two-loop functions 
with two external legs, including their infrared divergent on-shell derivatives. Results for one-loop multi-leg 
diagrams have been shown in [3j and additional material can be found in|3|. Two- loop three-point functions 
for infrared convergent configurations have been considered in |3j (hereafter III), two- loop tensor integrals 
in UJ. 

Many mass scales appear in the computation of physical observables within the Standard Model, gen- 
erating serious difficulties for the familiar analytical approach. Our purpose is to overcome these problems 
through a numerical approach. The application of our techniques has recently contributed to the evaluation 
of the two-loop fermionic correction to the effective electroweak mixing angle and of the full Higgs-mass 
dependence of the bosonic ones [7j. 

The approach described in J] is primarily intended for evaluation of multi-loop diagrams with internal 
massive lines. However, QED and QCD are integral part of any realistic calculation and they lead to infrared 
singularities. Therefore, any method aimed to a numerical evaluation of diagrams must be able to handle the 
infrared problem and infrared/collinear configurations should be treatable within the same class of algorithms 
used for the non-infrared cases or within some simple extension of the latter. 

For one-loop diagrams we have seen that our methods allow us to extract the infrared pole in dimensional 
regularization with a residue and a finite part that can be treated numerically [3]. The procedure has been 
extended in II to cover the on-shell derivative of two-point functions which are needed in the treatment of 
external legs. 

It is the purpose of this paper to extend the study of infrared divergencies to two-loop three-point 
functions. All diagrams are computed within the scheme of dimensional regularization [H] with space-time 
dimensionality n — 4 — e. Each loop in a diagram contributes at most one soft (zero gauge-boson mass) and 
one collinear (for zero fermion mass) 1/e term but the global order of the pole at e = can be greater than two 
due to simultaneous occurrence of ultraviolet poles which are removed by the introduction of counter-terms. 

To accomplish our goals we need an automatized procedure for handling infrared (and collinear) con- 
figurations: Landau equations [Hj represent the proper tool since a necessary condition for the presence of 
infrared divergencies is that the Landau equations arc fulfilled. Therefore, for each topology we build indi- 
vidual diagrams by filling all the lines with the line content of the theory, disregarding those configurations 
with vertex content not allowed by the theory itself. The generated result is examined and Landau equations 
studied for those diagrams that contain massless gauge-boson: if they are fulfilled then we have an infrared 
divergent configuration. The residue of the infrared pole(s) and the corresponding infrared finite part are 
then computed numerically. 

This part of the procedure is relatively easy while the difficult task is connected to the numerical evaluation 
of residues and of finite parts. They will be given in terms of multi-dimensional integrals over Feynman 
parameters with integrands that are not positive defined and, according to our strategy, their evaluation 
requires introduction of smoothness algorithms. 

Smoothness requires that, after suitable manipulations, the kernel in the integral representation and its 
first N derivatives be continuous functions and, ideally, N should be as large as possible. However, in most 
of the cases we will be satisfied with absolute convergence, e.g. logarithmic singularities of the kernel. This 
is particularly true when the large number of terms required by obtaining continuous derivatives of higher 
order leads to large numerical cancellations. 

There is a general approach for extracting infrared poles which goes under the name of sector decom- 
position We have examined this technique which, despite its great intrinsic possibilities, has its own 
problems: to name one it has been applied (so far) mainly to unphysical kinematics where infrared residues 
and finite parts are given in terms of positive definite integrands, i.e. it will not work properly around 
thresholds where the Feynman integrands are known to change their sign and imaginary parts show up. For 
recent developments see, however, ref. 

In our experience the form of the integrand, after many iterations of the sector decomposition technique, 
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is such that one can hardly imagine to design adequate smoothness algorithms. For this reason we have, 
quite often, privileged algorithms that keep under control the smoothness of the Feynman integrand at each 
step of the extraction of the infrared singularities. Usage of the whole machinery of hypergeometric functions 
has shown particularly useful in this respect. 

One may wonder why to devote additional efforts to the problem of computing infrared divergent dia- 
grams, given the spectacular success of analytical evaluation in QED/QCD: here we refer, in particular to 
the results by [12], by Q3| and by [Hj but also to [32]. 

The actual reason for pursuing this line of research is that QED and QCD are embedded in a more general 
theory, e.g. the standard model of fundamental interactions; from this point of view their handling is much 
more complicated. For instance there will be more than one mass scale for infrared divergent configurations, 
like in the decay of charged gauge-bosons and, with few exceptions, the analytical approach works only for 
very few scales or in the approximation where the scales themselves are arranged according to some fixed 
hierarchy, m <§C M etc. 

We are not claiming that a purely numerical approach is the final solution, rather one should carefully 
mix (semi) analytical extraction of dominant corrections (e.g. leading and sub-leading collinear logarithms) 
with numerical evaluation of sub-dominant, process-dependent, terms; the latter should be transformed in 
a way that allows for a safe, stable, integration where apparent singularities of the integrand are absent or 
limited to a minimum amount. Whenever a cancellation of dominant terms is foreseen we have to organize 
the calculation in such a way that these terms drop out before any numerical integration is attempted. From 
this point of view the technique of reduction of an arbitrary diagram to generalized scalar integrals is not 
always the best choice; master integrals quite often are individually more divergent (e.g. in the collinear 
limit) than the complete answer. Our technique does not grant any privilege to master integrals - from a 
computational point of view - and, therefore, seems more appropriate in handling the problem. 

The outline of the paper will be as follows: in Section [3 we define our conventions. In Section [21 wc 
review some of the tools that have been introduced to study infrared divergencies in quantum field theory. 
The connection between infrared divergent configurations and Landau equations is described in Section 13. II 
in Section 13.21 and 13.31 we present the procedure of sector decomposition while an alternative technique for 
extracting ultraviolet (if any) and infrared poles, based on properties of the hypergeometric function is given 
in Section 13.51 In Section 13.61 we present a discussion of threshold singularities. Starting with Sectional we 
present our results for all configurations, from Section lSTl to Section lSTl In Section l5~§l we present an explicit 
example of our procedure for classifying infrared divergent diagrams. Numerical results are summarized in 
Section [7] In Appendix we give a collection of technical details. 

2 Notations and Conventions 

Our conventions for dealing with arbitrary two- loop diagrams have been introduced in Sect. 2 of II. 
Conventions that are specific for three-point functions have been introduced in Sect. 2 of III; also the various 
families of two-loop vertex diagrams have been classified in III but, for the reader's convenience, they are 
repeated in Figs. I3HTT1 

2.1 Integrals and integration measures 

In particular, to keep our results as compact as possible, we introduce the following notations where 
xo = ya = 1: 




/ 



dC n ({x})f( Xl 



J ° 4=1 



(1) 
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Also, the so-called '+'-distribution will be extensively used, e.g. 



<IC,M:}) f 1 dx f{x,{z}) 
' x 



= / dC n ({z}) 



dx 



f(x,{z})-f(0,{z}) 



dC n ({z}) 



dC n ({z}) / dx 



dx 



x-1 
f(x, {z}) In" x 



X — 1 

f(x,{z})-f(0,{z}) 



= / dC n ({z}) 



dx 



In™ x 



The last relation in Eq.Q is used for evaluating integrals of the following type: 



dx 



dx 



f(x) Inx 



+ 0{e 2 ). 



(2) 



(3) 



Since we will have to split integrals during the evaluation of diagrams we have introduced a special notation: 



c,b 



dxF 1(S2 - a{x) = I dxFi- A (x 



dxF2- A (x). 



(4) 



In other cases we have to integrate over a triangle, for which we introduce the special notation (X = 1 — X) 



dx2dxz = 



(0,X ,xi) 



Xxi 



dx 2 



x 2 /X 



dx 3 + 



dxo 



X x x 



(xx—x 2 )/X 



dx 3 



(5) 



2.2 Alphameric classification of Feynman diagrams 

In our conventions any scalar two-loop diagram is identified by a capital letter (S, V etc, for self-energies, 
vertices etc) indicating the number of external legs and by a triplet of numbers (a, j3 and 7) giving the 
number of internal lines (in the qi , 52 and q\ — q2 loops respectively) . There is a compact way of representing 
this triplet: assume that 7 7^ 0, i.e. that we are dealing with non-factorizable diagrams, then we introduce 
K = 7max awx (/? — 1) + as — 1 +7 for each diagram. Furthermore, we can associate a letter of the alphabet 
to each k: for G — V we have a ma x = 2 and 7 max = 2, therefore, the following correspondence holds: 

121 -> E, 131 -> I, 141 -> M, 221 G, 231 K, 222 H. (6) 

For G = S we have a max = 2 and 7 max = 1, therefore 

111 -> A, 121 -> C, 131 J5, 221 D. (7) 

This classification is extensively used throughout the paper where we omit the suffix for scalar diagrams. 

2.3 Basic quadratic forms 

An x-dependent mass is always defined as 



x 1 — x 

In the following we introduce some quadratic forms that are widely used throughout the paper. 



X (x ; P 2 ; m , M) = -P l x 2 + (P 2 - m 2 + M ) x + m 2 = -P 2 [x - X % ) 
where we define 



B 



B^ = 



1 



X 4p2 



A(-P 2 ,m 2 ,Af 2 ), X 



P 2 - m 2 + M 2 — 



2P 2 



X x -1- x x . 



(8) 



(9) 



(10) 
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Here A(x, y, z) = x 2 + y 2 + z 2 — 2 (xy + xz + yz) is the usual Kallen lambda - function. If no ambiguity will 
arise we will simply write x(x). Furthermore, we introduce 

X(x; P 2 ; m, M) = x(l-x; P 2 ; to, M) = -P 2 x 2 + (P 2 — M 2 + to 2 ) a; + M 2 = -P 2 (x-X x ) 2 + B x , (11) 

/?(x, y 5 P 2 ; to , A/) = -P 2 x 2 + (P 2 - m 2 + M 2 ) xy + m 2 y 2 = -P 2 (x - X x yf + B x y 2 , (12) 

/3(x, y; P 2 ;m, M) = -P 2 x 2 + (P 2 - M 2 + to 2 ) x y + M 2 y 2 = -P 2 (x - X~ x y) 2 + B x y 2 . (13) 

Once again, we will drop irrelevant arguments if no ambiguity may arise. These quadratic forms will play a 
major role in the evaluation of the diagrams considered in this paper. 

3 Tools for virtual infrared divergencies 

Before starting a comprehensive study of two-loop, infrared divergent, vertices we collect in this Section 
a set of tools which are relevant for the general analysis of infrared divergencies in a spontaneously broken 
quantum field theory. 

First, we recall the classification of infrared configurations based on the study of Landau equations; we 
then move to a short review of the techniques which go under the general name of sector decomposition; 
the use of these techniques for handling integrable singularities in the evaluation of Feynman diagrams is 
also discussed. Alternative techniques based on a representation of infrared divergent diagrams through 
hypergeometric functions will be illustrated by means of simple examples. Mellin - Barnes techniques will be 
introduced to study infrared configurations around their normal thresholds and to extract collinear limits. 
Finally, we will present some of our new techniques which represent an extension of the Bernstein - Sato 
functional relations 17 . 

3.1 Infrared singularities and Landau equations 

As explained in the Introduction, the classification of infrared divergent Feynman diagrams is most 
conveniently based on the use of Landau equations P|. The whole procedure is better illustrated in terms of 
a simple scalar one-loop triangle. Consider the following Feynman parametric representation of the diagram: 




Figure 1: The scalar one-loop, three-point Green function with one massless internal line. All momenta are flowing 
inwards. 

^ Z" 1 f X l 2 

Co (pl,pl,P 2 ;m,0,m) = — dx dyVc 1 ~ e/2 (x,y), 

Z7r Jo Jo 

V a (x,y) = -pi x 2 - pjy 2 - 2px ■ p 2 x y + {pi - m 2 )x- (p% - P 2 - m 2 )y + m 2 . (14) 

To regularize infrared divergencies we use the scheme of dimensional regularization, with the space-time 
dimensionality n = 4 — e and where ji is the 't Hooft unit of mass. In principle one should carefully 
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distinguish between ultraviolet and infrared regulators (e > or e < 0) but for a comprehensive discussion 
we refer to |18| . 

It is important to recall that a necessary condition for the presence of infrared divergencies is that the 
Landau equations are fulfilled. A proper solution of the set of Landau equations for the one-loop vertex of 
Fig. n requires that 



P 2 = - 



If r I 1 / 2 ! 

— \2ml(m\ + ml)-{m\+ml+p\){ml+ml+pl)± X(-p 1 ,m 1 ,m 2 )X(-p 2 ,m 3> m 2 ) \, (15) 

which, for ni2 = and mi = = m is satisfied if p\ = p\ — —to 2 and P 2 is unconstrained . The 
corresponding V c reads 

V c 0, y) = m 2 (1 - x) 2 + to 2 y 2 + (P 2 + 2 to 2 ) (1 - x) y, (16) 

showing that V c = for x = 1 and y = 0. The parametrization in Eq. l|14fl is certainly not the most 
convenient one to evaluate the residue of infrared pole but this choice has been made deliberately, having in 
mind the more complex two-loop cases. 

In any software package designed to compute (automatically) physical observables the classification of 
the infrared cases must be a built-in procedure. The study of the Landau equations for a given family of 
diagrams is the most elegant way to classify its infrared divergent configurations. Again we refer to the 
one-loop vertex for illustration. The Landau equations are 

ai(g 2 + m 2 ) = 0, a 2 ((g + pi) 2 + TO 2 ) = 0, a 3 {{q + Pl + p 2 ) 2 + to 2 ) = 0, (17) 

Oil <7 M + a 2 (q+ pi)^ + a 3 (q+pi+ p 2 )^ = 0. (18) 

A solution of the system with en ^ 0, Vi gives the leading singularity. 

Let us multiply Eq. I|18|l by q, by p\ and by p 2 to obtain an homogeneous system of three equations where 
we use q 2 = —m 2 etc, from Ea. (|17|) . A necessary and sufficient condition to have a proper solution, i.e. not 
all the «i = 0, requires the determinant of coefficients to be zero, thus fixing a relation between internal 
and external masses. Any configuration that satisfy this constraint is a Landau singularity for the diagram 
which, however, does not necessarily imply that the diagram itself diverges at that configuration. For the 
Co function this corresponds to the well-known anomalous threshold. 

Consider the following case: p\ 2 = —m 2 ,mi i 3 = m and m 2 = 0. It is easily found that the anomalous 
threshold condition is satisfied, therefore the configuration is a Landau singularity; however, the question is 
of which kind. Let us insert the above values into the homogeneous system, what we obtain is 

™ 2 "i + {\p 2 + m 2 )a 2 = 0, (ip 2 +to 2 )c*i +m 2 a 3 = 0. (19) 

First of all we observe that it is not cither ct\ = or q 2 = — m 2 , etc; it can be both. Secondly, our 
configuration, where P 2 = (pi +p 2 ) 2 is unconstrained, is a singularity. Finally these is a special case of the 
general configuration discussed (P 2 free) which is even more singular giving, in the infrared case, the true 
leading Landau singularity. To have 0:1,3 ^ we must require P 2 = —4 m 2 which gives, in the annihilation 
channel, the well-known threshold singularity on top of the infrared one. This condition emerges also from 
the following argument: inserting p\ 2 — —m 2 in the anomalous threshold condition one obtains m 2 P 2 = 0, 
and to 2 (to 2 — P 2 — 4m 2 ) = 0, corresponding to the two signs of Eq. i|15|) . In a certain sense the constraint 
p 2 = — 4 to 2 is buried inside the anomalous threshold condition. To summarize: that all the propagators in 
a diagram are on-shell and that the consistency relation is satisfied does not necessarily imply that all di are 
different from zero: the infrared case is a clear example. Note the presence of a potential singularity also at 
P 2 = which, however, is not physical, i.e. is not lying on the physical Riemann sheet. The latter fact can 
be seen by inspecting the explicit analytical result, 

gQ R = ^ 1* |±jl+ IR.finite, ^ = 1 + ^- ( 2 «) 



our metric spacelike p implies positive p 2 . Further p4 = ipo with po real for a physical four- momentum. 



G 



Our strategy for the general classification of infrared divergent configurations will be on a diagram-by- 
diagram basis; we assume a certain number of zero internal masses with at least one unconstrained external 
momentum. Then we fix the remaining parameters to satisfy the consistency relation for the Landau equa- 
tions. Finally, we return to the original set of Landau equations and look for additional constraints that 
are necessary in order to generate the true leading singularity. The presence of a threshold-like singularity 
on top of the infrared poles is the sign that, after extracting these poles, we still have complications for 
the residual integrations that cannot be solved with a naive use of the method to be described in the next 
subsection. 



3.2 Sector decomposition 



Once a diagram has been identified as infrared divergent we must proceed with the extraction of the 
infrared polc(s) and with the evaluation of the infrared residue(s) and of the finite part. 

To continue our discussion we reconsider Ea. 1)14(1 and change variables, x — 1 — x' and y = (1 — x') y' . 
Next we split the integration by using 1 = 8(x — y) + 6(y — x), remap the two integrals into [0, l] 2 and obtain 



r~i 1 R 



dx 



dy 




-*- e )x,- W2 + 


(y~ 


1-e _ 


-xy-^Xy 1 ^' 2 


Xx - 


= m 


(1 


f yf + xy(xy~ 


2- 




+ P 2 y(l-x), 


Xy = 


- m 


"(1 


-f x) 2 + xy (xy — 


2 - 


2x) 


+ P 2 x(l - xy) 



The infrared pole can now be extracted, giving Cq r = P IR 

-l-e/2 



P -2-V- 

r IR — ^ . n 

I 7T 



dx dy x 



-P 2 y 



2 11 



= --—t d y 



m 2 (1- 



-yf 



P 2 y 



-l-e/2 



(21) 



(22) 



This procedure of extracting infrared poles, known as sector decomposition jl(J|. cannot be applied in its 
naive version to thee-point functions. The reason is that the residue of the infrared pole(s) is now given by 
an integral where the integrand is not positive defined over the parametric hyper-cube and, therefore, we 
cannot perform a straightforward numerical integration, at least if we want to avoid brute force and time 
consuming methods. Note that in our procedure we will use the van der Bij - Veltman parametrization 1!) 
of two- loop integrals and not the more familiar Cvitanovic - Kinoshita one [201 • 

To summarize the complete algorithm of sector decomposition, we perform the following steps (for each 
family of diagrams): At first, infrared and eventually collinear divergent configurations are selected by using 
the corresponding set of Landau equations which have been derived, diagram-by-diagram in III. Examples 
for two- loop vertices are given in Figs. 0^21 For a similar classification in QED with a massive regulator we 
refer to [21] • sector decomposition is applied to a graph G leading to an expression 



G = 



j2 s"( e ) 

1=0 



L 

E 

1=0 



Gi = 



1 



dxGi(x), 



(23) 



since each loop contributes at most one soft and collinear 1/e 2 term (e being the dimensional regulator) 
the highest value for L in a two-loop diagram is four; higher values in the Laurent expansion are only 
possible if ultraviolet divergencies are also present. In those cases the proper procedure is: ultraviolet 
poles are removed by local counter-terms, analytical continuation in e is performed and infrared poles are 
extracted. Smoothness algorithms are derived for each component G;, as indicated in Ea. 1(23(1 . where x 
is a vector of Feynman parameters, S is some simplex, Q\ is an integrable function (in the limit 8 — > 0) 
and Bf is a function of masses and external momenta whose zeros correspond to true singularities of G, if 
any. Smoothness requires that the kernel in Eg. ((23(1 and its first N derivatives be continuous functions and, 
ideally N should be as large as possible. However, in most of the cases we will be satisfied with absolute 
convergence, e.g. logarithmic singularities of the kernel. This is particularly true around the zeros of Bf 
where the large number of terms required by obtaining continuous derivatives of higher order leads to large 
numerical cancellations. 
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3.3 Integrable singularities and sector decomposition 



As it will become evident in the following sections, where we explicitly evaluate the diagrams, most of 
our integral representations will have the form 



G 



dx 77 \ 



In 



1 



A(x) 



B(x) 



or 



G(x) 



dx 77 \ 
s A (x) 



Li, 



fA(o 



\B(x) 



(24) 



where a; is a vector of Feynman parameters, S is some simplex and A, B are multivariate quadratic forms. 
Eq.J23Jl generalizes the Nielsen- Goncharov family of polylogarithms ^2S1 based on monomials in one variable, 
see Appendix Our integral representations are well-behaved around A(x) ~ but numerical instabilities 
could arise when simultaneously A(x) ss B{x) « as it will always be the case for collinear singularities. 
A nice solution to this problem |21| is to adopt a sector decomposition to factorize their common zero. 
Eventually, for some special configuration of internal masses and external invariants, this procedure will 
describe the correct behavior around a genuine singularity of G, e.g. A and B having a common zero of the 
same order. An example is provided by the following integral, 



J(a) 



f 1 1 


x 




/ dxdy — In 


1 + 


= J dxdy 


/o x 


L x + ay- 





In 



ay 



-In 1+ 

x V x + a 



(25) 



After performing a sector decomposition in Ea. l|25|) we obtain that a = is indeed an end-point singularity 
of J. All our numerical results are based on additional sector decompositions of integrands of the form shown 
in Ea. (|24() . although this is not explicitly indicated in the text. 



3.4 Infrared power-counting 

Before continuing the discussion on the evaluation of infrared divergent configurations it is convenient 
to introduce the concept of infrared power-counting. Given any parametric representation of a Feynman 
diagram we perform transformation of variables such that < Xi < 1 with 

/n n 
dC NX -/, X N = E a h-in II x i 6 (E l i - L ) ( 26 ) 

h...l n i=l i=l 

and Xn = 0, for xi — 0, Vi, with L > 1. After introducing polar coordinates the integral can be written as 

G = J dfi n _! ^ drr^^x^^n), (27) 

with Xw,rod(0, Vt) ^ 0. Define A, B such that n — [iL— 1 = A — Be. IfA>0 the integration does not lead to 
an e-pole and the corresponding diagram is infrared safe. A typical example is represented by V E of Fig. [3] 
where, although x E = at the hedge of the integration region, power-counting shows that no infrared pole 
arises. 



3.5 Infrared singularities and hypergeometric functions 

We have seen in the previous sections that a general method for extracting the infrared poles (within 
dimensional regularization) of Feynman diagrams is based on sector decomposition. However, in most of 
the two-loop cases sector decomposition has drawbacks: first of all the number of generated sectors tends to 
increase considerably and then the procedure creates new integrands with polynomials of very high degree, 
the higher the number of iterations the higher the degree. The consequence is that one cannot find an 
adequate smoothness algorithm to handle the final integration. We give a second, alternative, procedure 
once again illustrated in terms of a Co function. 

The polynomial V c of Ea. (|14|l can be rewritten (P 2 = — s, pf = — m 2 ) as 

V c (x, y) — m 2 x 2 + m 2 y 2 + (s — 2 m 2 ) x y — 2 m 2 x — (s — 2 m 2 ) y + m 2 , (28) 
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and we perform the transformation y = y' + ax with to 2 a 2 + (s — 2 to 2 ) a + m 2 = 0. This transformation 
is designed to make V linear in x and, whenever possible we always seek for a transformation that makes 
the Feynman integrand linear in one of the variables. After some straightforward manipulation we obtain 
i n 2 C = /i e (C 1 + C 2 ) weher 



C£ = (l-a) J dy I dxV- 1 t/2 (x,y), C 2 = a / dy 
with ^-polynomials given by 



dx v 9 



-l-e/2, 



(29) 



= 
V 2 = 



s (1 + a) y — s + 2 (1 — a) m 2 i-asy 2 -2(l-a) m 2 y + to 2 , 



a s + 2 (1 — a) to z 



«s-(2a-l) m 2 



In Ea. 129fl the x-integration has the general form 



h(y) 



dx 



B i {y)-A i {y)x 



-l-e/2 



2^(1 + 6/2, 1; 2; ^y). 



1,2. 



(30) 



(31) 



Using well-known properties of the hypergeometric function we obtain, for | arg(— z)\ < ir, 



a Ji(H 2 . . . _. 



1; 2; z) 
from which we derive 



2 Fx(l, 1 + -; 1 + -; l-z)+ (1 - z)^' 2 2 F X {\ , 1 - - ; 1 - - ; 1 - z)\ (32) 



e/2 



-e/2 



Given the form of in Ea. H30[l we can simply expand around e = obtaining 



C x Q ={a-l) / dy 



My) 



in 



Ai(y) 



which is well-behaved for A\(y) = 0. However, for i — 2, we find 

A 2 (y) = a(s,m 2 )y, B 2 (y) = b{s,m 2 ) y 2 , 



(33) 



(34) 



(35) 



a(s,m 2 ) = -as + 2(a-l)m 2 = - (0+ 1) 2 to 2 , 6(s,to 2 ) = -as+(2a-l)m 2 = ^(/3+1) 4 to 2 . (36) 



Therefore, C 2 is infrared divergent and we get 



Cl = -^b^' 2 { S ,m 2 ) / dyy" 1 - In 



a(s, to 2 ) 



= 2 



a ,0-1 
In 



a(s, to 2 ) /3 + 1 



6(s, TO 2 ) 
1 1 1 P~ l 1 1 U 2 \ 



1 - - In 

4 



a s, to 



6(s, TO 2 ) 



(37) 



with f3 2 = 1 — 4m 2 /s. The infrared pole has been isolated and infrared residue and finite part are already 
in a form that allows for direct numerical integration. For a general diagram they are rational functions of 
the residual Feynman parameters and we have been able to derive adequate smoothness algorithms for their 
integration. 



3.6 Threshold behavior and Mellin-Barnes transforms 

In computing infrared residues and finite parts of diagrams we face one additional complication: in 
general we end up with integrands which are not positive definite. Therefore, around those configurations of 
the external parameters where threshold singularities occur the algorithm has to be modified. The general 
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idea is to isolate the singular behavior and to write, for the regular part, an expansion in some Kallen 
function of the external parameters. We borrow the relevant technique from another problem, the large 
energy expansion of Feynman diagrams which is best performed by using Mellin-Barnes transforms |22| . 
The whole idea is better illustrated with a simple example. Consider the infrared residue of Cq: 

P IR = --J^L B 1+e/2 (-m 2 ; to 2 , P 2 + Am 2 ) , 
1 



B a (p 2 ; m\ , to 2 ,) — i dx x( x ) — id , x( x ) = — P 2 x 2 + (p 2 + m \ — m i) x + m ii (38) 
Jo L - 1 

and <5 — > + . Let a;± be the roots of \ — i $ = 0. If they are complex or real but external [0, 1]) then the 
numerical evaluation is straightforward; when they are real and internal (6 [0, 1]) the integration contour 
can be distorted and the integral can be computed unless a pinch will occur in [0, 1] |23| (or an end-point 
singularity); this happens for A(— P 2 , to 2 , to 2 ) = 0. 

The integral in Eq. (|38() is simple enough to be computed analytically but in our approach we will pretend 
to treat it numerically; distortion is performed, unless a pinch occurs. Then the question will arise of what 
to do around those parametric regions where A = 0. Here we describe our solution for extracting the leading 
and sub-leading behavior of P IR around A = 0. First we rewrite x as 



x(x) = to 2 (x - x ) 2 - - iS, x = -^-l, (39) 



4 TO 2 " 2 TO 

and consider the case P 2 — — s with s > 0. Further, we assume that < xq < 1 and obtain 

-2 a pi r+ioo \ 

B a = dx ds B(s,a-s)p a - s Q- s , p~ x = r-iS, Q = (x - x ) 2 - i 8. (40) 

2m ,/ J_ ioo 4to 4 

where B denotes the Euler beta- function. Ea. (|40|l is valid in the vertical strip < Res < a. We choose 
a < 1/2, require < Res < 1/2 and perform analytical continuation to obtain 

1 1-150 yl-2s 

d *Q- s = E y^Ys- (41) 

X=xo 







Since we are interested in the limit \p\ — > oo, the s-integral will be closed over the right-hand complex 
half-plane at infinity, with simple poles at s = 1/2 and s = a + kk>0. In this way we obtain 

11 / A \l/2— a 

B a = B(~,a- -)mr 2a + sub-leading, s^4m 2 . (42) 



'2' 2' V 4to 

The result can be easily generalized to the case of unequal masses. In this case we have also have a pseudo- 
threshold, s = {mi — TO2) 2 , which is not a singularity on the first Riemann shhet since Xq $ [0, 1]. 



3.7 Old and new Bernstein - Sato - Tkachov functional relations 

For all one-loop multi-leg diagrams we have developed computational techniques based on the proposal 
introduced in I. Most of two-loop infrared convergent diagrams can be computed by following the same strat- 
egy. However, some extension has to be introduced to deal with the general case. To derive new algorithms, 
we recall the definition of Bernstein - Sato polynomials |17j : if V{x) is a polynomial in several variables then 
there is a non-zero polynomial b(p) and a differential operator V(fi) with polynomial coefficients such that 

V(jjl) V» +1 {x) = b{p) V»{x). (43) 

The Bernstein-Sato polynomial is the monic polynomial of smallest degree amongst such 6(/i). If V(x) is a 
non- negative polynomial then V ll {x) 1 initially defined for fi with non-negative real part, can be analytically 
continued to a meromorphic distribution-valued function of p, by repeatedly using the functional equation 

V"(x) = -^V(n)V* +1 (x). (44) 
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The Bernstein - Sato - Tkachov theorem [25] tells us that for any finite set of polynomials Vi(x), where 
x = (xi, . . . , x N ) is a vector of Feynman parameters, there exists an identity of the following form (hereafter 
a BST identity): 

V (x,d)]J Vt + \x) = B v [] Vt(x). (45) 

i i 

where V is a polynomial of x and di — d/dxf, B v and all coefficients of V are polynomials of ^ and of the 
coefficients of Vi(x). Furthermore, if the polynomial V is of second degree we have a master formula, due to 
F. V. Tkachov [23]. We write the polynomial as V{x) — x l H x + 2 K f x + L, where x t — (x\, ...,x n ), H is an 
n x n matrix, K is an n vector. The solution to the problem of determining the polynomial V is as follows: 



V = 1 



(x-X v ) t d x 
2(^ + 1) : 



B v = L- K l H^ 1 K, 



H^ 1 K. 



Therefore we have: 



1 

~B~v 



(x-x v ya x 



2(M + 1) 



V^ +1 {x) : V~ l {x) 



B v 



l--{x-X v ) t d x \nV{x) 



(46) 



(47) 



The list of BST relations must be extended to cover infrared singular cases. It often happens that V(x) is 
not complete and so other BST relations, originally defined in Ea. H45|) . are needed. A typical example is 
when V is linear in one variable V(x, y) — h {x — xq) 2 + c y + b. For this polynomial the H matrix is singular 
and it can be easily shown that the following relations hold: 



= - 
b 



1 



(x - x )d x + 2yd y 
2(/x + l) 



V" = --^d y V^ +1 . 



(48) 



c /i + 1 

A better way to proceed without loss of generality is to introduce a Vo and a V\ with the property that 

V±=V Q ±V\d x , V+V(x)=B v . (49) 

BST relations can be written as 



1 
~B~v 



V 



1 



V\d x V» +1 (x), for 



V~\x) 



1 

~B~v 



V + V\ d x In V{x) . (50) 



A first extension of Eq. (|47|) is given by the following example 

1 



V-\x) ln n V(x) 



1 

~B~ V 



V + V\ d x In V(x) \n n V(x). (51) 

n + 1 

Eq.JHU can be easily generalised to arbitrary powers of V: 

W{x) In" V(x) = ±- {Po In" V(x) + £ 7-7^+1 V \ d * ^ V^\x). (52) 



As a next step we seek for similar relations, holding for Nielsen polylogarithms and for the hypergeometric 
function These new relations can be obtained with the help of the following formulae: 



Li„ + i (x) = - Li„ (a;) , 
dx x 



d a — 1 

2 Fi(a~ 1 , b; c; x) 



dx 

Given a generic polynomial A{x), it is easily verified that 



2 Fi(a , b; c; x) - 2 Fi(a- 1 , b; c; x) . (53) 



V 



■ Li„ 



1 

~bZ 



A 



LL 



V+ A - V\ d x Li n+1 



(54) 
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V- a 2 F 1 (a,b;c; = 1- 

n 



]^I[ 2 Fi(o,6;c; £>+ A - a *i(a - 1 , 6; c; ^)V{d x A 



1 - a 



^-^FiCa-l, 6; c; 



-4. 
V 



(55) 



In particular for n = 0, 1, we have Lio (x) = x/(l — x) and Lii (x) = — ln(l — x), for which the corresponding 
relations are 





1 A 


1 






V V - A 


~ 57 


V 


1 




1 


1 










V 




= 


A 



1 



In 1 - 



V + A + V{d x U 2 



When the BST factor, B v , vanishes we have: 

1 



V — A 



V + A = - V\ d x In 1 - 



for B v = 0. 



Another set of (new) relations is obtained starting from known properties of poly logarithms, 



d_ Li w+ i (x) 
dx x 

For x = A/V, when n = 1, we easily obtain: 



Li„ (x) - Li n+ i (x) 



y2 



In 1 - 



1 



Using Ea. (|60|l for n — 2, we obtain: 



.4 



V 2 



V 



A 



AV-- ln(l- 



V 



V\ d x In 1 - 



(56) 
(57) 

(58) 

(59) 

(60) 



775 Li 2 -T =^MO_-Li 2 - 



I? 



.4 



.4 



V 



V 



— In 1 - — 



.4 



1 -y 



(61) 



The same procedure can be iterated to get a similar relation for A/V 2 Li„ (V/A). 

One can derive other relations of this type but here we have restricted our attention to those which are 
actually used in this paper. Their application in computing Feynman integrals is based on the integration 
by parts of the terms containing V\ d x (when it is not applied just on ^4). 

After applying BST relations we end up with integrands which show a less divergent behaviour with 
respect to the original one; in most cases the integrand is smooth enough to allow for a stable numerical 
integration. When this is not the case we reiterate the use of BST functional relations; they are followed by 
a second integration by parts, etc, etc. 

A typical example of smoothness algorithm corresponds to apply a BST functional relation to Eq. l|22|l . 
In this case we obtain 



i 7r 2 m 2 P IR 



(EY / d yv- l -'/\y), v(y) = y 2 + (z + 2)y + l, z 

Vm/ In 



2 



p2 



(62) 



dy v 



-l-e/2 



(y) 



1 - - (z + 4) In (z + 4) + - 
4 2 



dy In v(y) + O (e) 



(63) 



showing the additional threshold singularity at P 2 = —Am 2 . Note the absence of a singularity at P 2 = on 
the first Riemann sheet. 
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4 Derivatives of two-loop self-energies and infrared poles. 



In II wc have defined the on-shell derivative of a two-point function, where possibly some of the internal 
masses are zero, as the d/dp 2 derivative evaluated at the mass shell of one of the non-zero internal masses. 
These derivatives are used to construct wave function renormalization factors and a complete list of results 
has been shown in II. 

Here, we briefly review the subject; consider the on-shell derivative of S° (for the simplest topology 
S A the on-shell derivative is infrared finite). Once again, a necessary condition for the presence of infrared 
divergencies is that the Landau equations are fulfilled: for S c , we see that s = (mi + m% ± m^) 2 and 
m\ = {mi + rn2) 2 are satisfied by rri2 — TO4 = 0, mi — TO3 — m and s — m 2 . It is easily seen that S° is not 
infrared divergent but its derivative with respect to p 2 shows an infrared pole when computed on-shell. To 
continue our discussion we consider the case mi = TO4 = and mi — 7713 = m, a typical example of which 
is shown in Fig. After a straightforward calculation we obtain 




Figure 2: A two-loop diagram contribution to the VK-boson self-energy. 



^ (-£»)' 

m \ 7T m z I 



2 1 

4-2(2-7)- 



7 + 7 (4_ 7 )_I C (2) + 0( e ) 



(64) 



Results for S® = S^ 31 have been derived in Sect 7.4 of II; results for S° = S 221 have been derived in Sect 
7.9 of II; 7 is the Euler constant and £ the Riemann zeta function. 



5 Infrared divergent two-loop vertices 

Starting with this section we present a detailed discussion, diagram - by - diagram, of two-loop infrared 
divergent vertices. Once again, evaluating a specific diagram means to derive a set of algebraic manipulations 
that return a multi-dimensional integral representation of the diagram for which numerical methods can be 
safely applied. For infrared configurations this means evaluating both the residues of the infrared poles and 
the finte parts. Special limits always require additional refinements of the procedure. The key ingredients in 
our derivation are: BST functional relations, the whole machinery of hypergeometric functions and Mcllin- 
Barnes transforms. 

Several diagrams belong to the general G lNl family which is specified in terms of a set of momenta fc, 
which are linear combinations of the external momenta pj ; fc^ = pi + . . . + p^. Our parametrization is as 
follows: 

G lNl = -(^j T(N-2 + e)j\x J dS N (y,u 1 ,---,u N „ 1 )[x(l~x)Y* /2 (l-yY^xl-r^ 

XlN1 = u l Uu + 2K}u + {m 2 x - m 2 3 ) (1 - y) + m 2 , (65) 
where we have introduced the following quantities: 



n, s /',-/',• >Ci = ^ (kf - kf_i + m 2 i+3 - m 2 i+2 ), i,j = l,...,N—l (66) 



and where the ^-dependent mass is defined in Eq.©. Three out of six of the non-trivial two-loop vertex 
families can be described in terms of the representation of Eq. I|65|) . 
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5.1 The V E diagram 

The V E diagram of Fig. [3] is representable as 



7T 4 V E 



d n qi d n q 2 



[1] E [2] E [3] E [4] E ' 



(67) 



[l) E = ql + ml, [2} E = (qi-q 2 ) 2 +m 2 2 , [3] E = [q 2 - p 2 f + m 2 3 , [4] E - (q 2 - P) 2 + mj. (68) 
Although V E does not exhibit infrared poles it is simple enough to illustrate our procedure. A necessary 




Figure 3: The irreducible two-loop vertex diagrams V E . External momenta are flowing inwards. 

condition for the presence of infrared singularities is that the corresponding Landau equations are fulfilled. 
For V E this requires the condition A P 4 + 2 B P 2 + C = 0, with coefficients 



.4 : 
C 



'3- 



4 777-3 



D 



p 2 + (mi + m 2 y 



(pi - m i + m l) - 2 m l (pi + pI), 



B 2 



X(-pi, ml, ml) A(-P2, ("7i + m 2 ) 2 ,ml) 



(69) 



Let us concentrate on infrared divergencies due to photons. Then any tri-linear vertex has at most one 
photon line and any quadridinear one has at most two photon lines. We are looking for solutions of the 
Landau equations where two external momenta are on some mass-shell and the third one is unconstrained. 
The procedure will be as follows: 

Let mi = and P 2 be a free parameter. We immediately find that Landau equations admit a solution with 
at 7^ 0, Vi if 7773 = 0, pf = — m|, p 2 = — m\. 

Let mi = and p\ be free, we find 777,2 = 0, P 2 = — m|, p 2 = — m|. 
Let mi = and p\ be free, we find 7774 = 0, P 2 = — m|, p\ = — m|. 

Another possibility is to start with 7713 = 0. We proceed as follows: 

Let m.3 = and P 2 be free. We obtain p\ = — m\, p\ = — (m% + m 2 ) 2 , which is of no interest since we do 
not expect a theory with such a peculiar relation among masses. 

Let 7713 = and p\ be free. We obtain mi + 77i 2 = 0, which is of no interest because it requires three photon 
lines in the same vertex. 

Let rri3 = and p 2 be free. We obtain 7774 = 0, p\ = 0, P 2 = —(mi + m 2 ) 2 , again of no interest. 

Because of the symmetry of the diagram the m 2 = and 7774 = cases are already covered. Consider 
mi = 7773 = and pf = — m 2 ,p 2 = — m 2 . Using Ea. (|65|) we obtain 



X E (x,y,z) = -ml 



i-y 

l-x 



- ml (1 - yf ~m\z 2 - (P 2 + m\ + m 2 ) (1 - y) z + m\ (1 - y), 



(70) 
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showing that x E = for the specified configuration at y = 1 and z — 0. However, x E appears with exponent 
— e in the parametric representation for V E so that no infrared pole will show up. To continue our analysis 
we consider the following cases: 



ii = — and P 2 free. We obtain p\ 



-ml, p 2 



mi = m 3 = and p\ or pi, free. We obtain that m 2 = and the case is of no interest because it requires 
three photon lines in a vertex. 

The last case is m\ = 7773 = 7774 = which requires either p\ = or 7772 = 0, i.e. more than two photon lines 
in a vertex. 



5.2 The V diagram 

The V family of diagrams, shown in Fig. is representable as 

1 



with propagators 



ir 4 V = p 



2< 



d n qi d n q 2 ■ 



1],[2],[3],[4],[5],' 



(71) 



[1], = q\ + ml [2], = (gi - q 2 f + m 2 2 , [3], = q\ + ml 



Wx = fe+Pi) 2 +mi [5}r = (q 2 +P) 2 +ml (72) 
Referring to Section l3~Tl we see that the leading Landau singularity for V is equivalent to the sub- leading 




Figure 4: The irreducible two-loop vertex diagrams V 1 . External momenta are flowing inwards. 

one, with a\ — a 2 = 0, when the peculiar condition 7713 = 7711 + 7772 is added. Therefore we discuss first 
those infrared configurations that correspond to sub-leading Landau singularities, and postpone the special 
case 777,3 = m i + m 2- 

If «3 = the reduced diagram, i.e. the one where the third propagator is shrunk to a point, corresponds 
to a V E topology which is free from infrared poles, as explicitly shown in Section I5~T1 

If a% = 012 = the reduced diagram is a one loop three-point function and the classification of the infrared 
singularities is simpler. We obtain: 1) 7773 = 0, P 2 



or 3) 7775 =0, P 2 = — m§, p\ = —m\ 
the third are actually the same). 



-to§, Pi — —7774, or 2) 7774 = 0, p\ 



22 2 
-777,3, P2 = ~ m *' 



They correspond to the configurations shown in Fig. (the first and 



5.2.1 Extraction of UV and IR poles 

For both configurations 1) = 3) and 2) we perform the transformation y = 1 — y'\ the diagram reads as 
follows: 



V/ = - 



r(l + e) f dx f dS 2 ({z}) f 1 dylxil-x^-^y^-^a^y + Z^y 
Jo J Jo 



(73) 
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P2 



Pi Pi 

Figure 5: The V 1 infrared configurations. The photon line represents a general massless particle while the dashed 
and the continuous lines represent different massive particles. The mass of the two particles in the bubble are mi 
and rri2. 



where i — a, 6, and a a;/ = m 2 , a^-j — m 2 — M 2 , with m 2 defined in Eq.JSJ), and 

Z a ;r = /? (z2 , zi ; P 2 ; to 2 , M 2 ) Z blI =p(l- Zl ,l-z 2 ;P 2 ;m 2 M 2 ) (74) 

and (3 defined in Ea. l|12[l . V is also ultraviolet divergent and we look first for a procedure that extracts the 
ultraviolet pole. Consider the innermost integral appearing in Eci. (|73[) 



y., • dyy^(l + ^y)-^. (75) 

It is convenient to evaluate this integral in terms of an hypergeometric function |27|: 

ViM = ~ Z-J- e 4 /2 2F l (l + e, e -;l+ e -;-±), = (76) 

e 2 2 Q-j ai-j{l-zi) 

Using a well-known property of hypergeometric functions (see appendix we obtain 

^ ~ ~ r(l + e) t;I 4;I ~ 2T^ i;I ( ^ 2^1(1 + e, 1 + - , 2 + -,-Ci;i)> (77) 

with the result that the ultraviolet pole has been extracted. The second term in Eq.JTTJ) is finite, so we can 
set e = 0, obtaining for V/: 

V i=-(^) e - e T2 { 1+e i) X * Z <+f Q dx J dS 2 ({z}) l ^±^l+0(e), (78) 

Xi- tl = J dx [x(l - x)]-c/ 2 a - e/2 , Z, = J dS 2 ({z}) Z^ t/2 . (79) 

In the last expression we have dropped out the index i for Z Iy since both cases produce the same result (it 
is easily seen by transforming z\ —* 1 — Zi, z<i — ► 1 — Z2 and z\ <-> z 2 , for the case b). 

The second term in Ea. l78|) is well-behaved when integrated over x, z% and z 2 and additional manip- 
ulations are not needed. As far as the first term is concerned we can write, after some straightforward 
manipulation, 

e r 1 e 2 r 1 

X i;I = l-~ / dx\nV i , I + - / dx\n 2 V V:I (80) 



2 ./,, ' ' n 

where the two new quadratic forms are 

,2 / 1 \ 1 2 



V a;I = m\ (l-x) + m\ x, V b;I = X (x ; —M 2 ; m 2 , m 2 ) . (81) 
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For Zj we map the integration region into the square [0, l] 2 obtaining: 

Z,= JdS 2 ({z})(3- 1 -^ 2 (z 2 ,z 1 ;P 2 ;m 2 ,M 2 ) = JdC 2 ({z})z^ X ~ 1 ^ 2 {z 2 ) = --J 1 dz X - 1 -^ 2 {z). (82) 

where (3 has been defined in Ea. (|12l> . For the ^-integral we use one BST iteration, integrate by parts and 
expand in e obtaining: 

Zj =TZ[e- 1 + 1Z' 2 +Tl I 3 e, (83) 
where the coefficients in the e expansion are 



2™ nl B, 



dz In" 1 x( z ) ln%(z) + 2n 



= X In" + X In" x(0). (84) 



The BST factor B x and co- factor X x are collected in Eci. (|10fl . Of course, the BST method fails when masses 
and external momenta are such that B x is very small. The general solution to this problem is based on 
the method of Mellin-Barnes transforms. For this particular diagram, however, the situation is very easy. 
Referring to Eci. (|82J) we have integrals of the form 



Hi, Hj 



dzxj 1 e/2 (z), xA z ) = x( z ; p 2 ; m l, m l) 



(85) 



which we want to evaluate in the limit B x — > (x is defined in Eq.JJjJ). The integral in Ea. (|85|) is a B a 
function of Eq. l|38[) and can be evaluated for an arbitrary values of the exponent a. Setting a = 1 + e/2 in 
Ea. p%|l we obtain 



Bi+e/2 (p 2 ; m 2 a , ml) 



B 



1 1 + e 



2' 2 



-,(1+0/2 



X=a 



-I-e/2 



X=a * Z L P X=a 



lnX 2 s + 2 
X 



o( P - 



where B is the Euler beta function and where we have introduced 



4s 2 



P 



s, 



s + mz — mt 



A = A(s, to 2 , to 2 ), 



(86) 



X-i6' 2s 

where X(x, y, z) is the Kallen function and where we assume s > 0. Collecting all pieces together, we get 



(87) 



V% = 2K[, V% = 2K I 2 -K I 1 [ dx In V i;I , 

Jo 



Vy =2K 



K - j^dx K* 2 InVi,, - ^ In 2 V v , - J dS 2 ({z}) 



ln(l + Ci;/) 



Zi-r 



(89) 



5.3 The V M diagram 

The V M family of diagrams, given in Fig. is representable as 

1 



7T 4 V M = /i 2e 



d n qid n q 2 



[1]m[2]m[3]m[4] m [5]m[6]m' 



(90) 



[1]m = q\ + m\, [2] M = (qi - q 2 ) 2 + m 2 ., [3] M = q\ + to 2 , 
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Figure 6: The irreducible two- loop vertex diagrams V M . External momenta are flowing inwards. 



[4] M = (q 2 +Pi) 2 +ml [5] M = (q 2 + P) 2 + m% [6] M = q\ + 

If 7TI3 ^ ttiq then is the difference of two V diagrams, 

1 



V = 



V'(P ;mi,m 2 ,r7i 3 ,m4,TO 5 ) - V'(P ; mi, m 2 , m 6 , m 4 , m 5 ) 



(91) 



(92) 



In general Landau equations are the same for V M and for V so that the classification of potentially infrared 
singular V M vertices follows closely the discussion presented in the previous section (see Fig. , the main 
difference being in the exponent of the integrand. Configuration c) is a special case of configuration b) and 
will be treated in a second step. 

5.3.1 Evaluation ofV M , cases a) and b) 

For the first two configurations of Fig. [7] we have 



r v 



1-21 



T(2 + e) dx dS 2 ({z}) / dy[x(l-x)}-^ 2 y^ 2 - 1 (l-z 1 -y)(a i , M y + Z i , M y 



2-e 



(93) 



where i — a, 6, a a - M = m 2 > ab-. M = rn 2 , — M , with m 2 defined in Eq.JB|) and with 

Z a -M = fife, z 1 ; P 2 ;m 2 , M 2 ) Z b , M = (l - z x , 1 - z 2 ; P 2 ; m 2 , M 2 ) (94) 

and defined in Ea. p2(l . V M is also ultraviolet divergent and we need, once more, a procedure for extracting 
the ultraviolet pole. Consider the innermost integral in Ea. (|93l) 



Vt l = - 



0£\ r(2 + e) J dx J dS 2 {{z})[x{l~x)]-^ 2 y^ 



yi;M — z. 



1- Zl 



d yy ^(l-z 1 -y)(l + ^y)- 2 ^. 



The integral of Ea. (|95|l is expressible in terms of hypergeometric functions, 



f b 1 

/ dyy s ~ 1 (l + ay)~ u = -b^F^v, s; s + 1; -ba) 7 Re s > 0, | arg(l + a b)\ < n. 
Jo s 



(95) 



(96) 



The result is further transformed according to well-known properties (see appendix Q and we obtain that 
1 



r(2 + e) 



7 -l-e/2 -c/2 
^i-.M a i;M 



(97) 
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Pi 



Figure 7: The V infrared configurations. The photon line represents a massless particle while the dashed and the 
solid lines represent different massive particles. Furthermore, the masses of the two particles in the bubble are m\ 
and rri2. 



where the function $j ;M is defined by the following relation 

2 



2 



(i-zir 1 ^ 2 



4 £ 2^(2 + 6, 2 + ~; 3 + |; - 



2J F 1 (2 + e,l + |;2+|;- C ^) 



1 



_ ln(l + Ci; M ) 
Cj;m 



(98) 



where we used Ci;m = Zi- M /ai- M (^ — zi). In Eq.JHHJl we have put e = and the last result follows after using 
the well-known relations for the hypergeometric function collected in appendix IU1 
After inserting these relations in Eq. l|95[l we obtain 



V M = - 

where we have introduced 



2 ^ e / f \ r/2 



(99) 



XI 



o 



dx[x(l-x)]~^a]-r e/ \ Z\u= dS 2 ({z})(l- Zl )Zr 



-2-e/2 



dS 2 {{z})Z\ 



-l-e/2 



(100) 



Note that X}. M has been already computed while analyzing the V diagram, Ea. (|80|) and Ea. (|81|l : therefore 

Xi- M = X i;I , V a ,b;u = V a ,b-,i, (101) 



19 



For X 2 M we distinguish between the two cases a) and b). The computation of X 2 . M is straightforward, since 
o-a;u is positive definite in a; € [0, 1], and gives 



x 1.m = x lt + X i-L ^ i = a , fo ; 



(-1)" f 1 



v2,n 



dxx{l-x)V-" ln"V a;1 



(102) 



For X 2 . M instead we use the BST-method to increase the power of 14 ;M and obtain: 



X, 



2.0 

b\M 



1 



273 A 
where we use 



X, 



2B A 



dx < In Vb-, 



x(l-x)--^lnVb;i 



(103) 



3x 2 -2{l + X M )x + X M , 



X M — 



M 2 +m\-ml 
2M 2 



Br 



A(A// 2 ,to 2 ,77tJ;) 
4M 2 



(104) 



The case B M = requires M — mi±m 2 which in all realistic situations means 7711(7712) = and 777,2(7711) = M. 
This is equivalent to configuration c) which will be treated in the next subsection. 

When examining Z\. M we transform z 2 — > Z\ z 2 for i = a and Z\ — > 1— Zi, z 2 — > 1 — 22, Zi Z2, #2 — ► 21 Z2 
for i = b. In both cases the z\ integration can be performed, giving 



4m = / dC 2 (l-z 1 )z 1 - 3 - e X (z 2 ) 



-2-e/2 _ 



(i + c)(2 + c) y 



4 M = / dc 2 zr 2 ^z 2 x(^)- 2 ^ /2 



1 



l + e 



dz 2 z 2 x(zz) 



dz 2 x{z 2 ) 



-2-e/2 



-2-e/2 



(105) 



For the ^-integration we use twice the BST-method followed by integration by parts and expand around 
e = obtaining: 

1 



4R 



/ T~t i ; M 1 % ; Af 



i = a,b 



(106) 



2R 



S^ 1 - X x In M - X x In m 2 + 2 



M 2 m 2 



1_ /InAf 



73, 



M 



v~x ■ —x • ~> • 8 - x 
X - ' -01 - - • - V 



-772 (lnM' + 4) + -X 



1 f\nm 2 4 \ 2 ,s 
2 (mm 2 +4) 



8 x V B x ,n 



K' M = - 75^ (*x " *x ^ - *x + 2) - 2 -I 



7^ ;M 



1 



4/3, 



X S*° 2 - 85"" + 14 X S* 1 - X 2 lnm 2 (lnm 2 +6) - X InM 2 (X InM 2 + 61-4) 



x 

4+12X, 



+ ^ (InAf 2 + 3), S™= f dz z n ln m X {z). 



(107) 



For .2^ we can follow the derivation given for V: 



zl = nr e^ 1 + n 



2 ■ 



7? M = u 1 



(108) 



with 72-^ given in Eq. l|84|l . One expression remains to be computed: 



Jr. 



dS 2 {{z}) / dx 



1 - 



ln(l + C< ;M ) 
Ci;*f 



(109) 



Here there is an additional problem, namely the integrand is well-behaved in the limit Zi- M — > but not 
when a,i-M —> 0. 
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For case a), a a - M never vanishes in the integration interval and we can simply write: 



Ja: 



dC 2 ({z}) / dx 



ml zix{z 2 ) 



ln(l + i] a . M ) 



TJa;M 



ZIX{Z2) 

ml{l-zi)' 



(110) 



On the contrary, to compute Jb-u is more convenient to reexamine Ea. (|98|l . Setting e = we have: 



J,, 



b:M 



dS 2 ({z}) I dxa^i I (l~z 1 )- L [- 2 F 1 (2,2;3 ] -C b ;M)-2F 1 (2 1 l;2;-C b , M )\ (111) 



As the next step we use the BST relation of Eq. l|55|) which, in the present case, reads: 



A 



A 



Vr^ Fl (2,b- c; p^) = V b -^ Fl {2, b; c; ^L) 



x - X, 



Vb-M B 

x X M Q x A\)- M 



d x 



A 



Vb-.M 



vr 1 



2 F 1 (2,6;c;^)- 2 ^ 1 (l,6;c;^)ll. (112) 

Vb-.M Vb-M J 



where B M and X M are defined in Ea. i|104[l and A^ M = — x (1 — x) Zb- M /(l — z\). After integration by parts 
we get: 



Jb:M 



1 



2 B M 



dS 2 ({z}) 



1 dx^—hx{l-x) \l - ln(1 + Ct;flf) 1 +a M ln(l + C i;M )j (H3) 



where a M = 3x 2 — 2 (1 + X M ) x + X M as defined in Ea. (|104l) . To obtain the last formula we have used the 
explicit expression for the hypergeometric functions given in appendix IU1 

The case where masses and momenta are such that B x is vanishing small must be treated separately. In 
the following we will show how Mellin-Barnes techniques can be used in this kinematical situation (for the 
cases under discussion). For case a) the situation is particularly simple: from Eci. (|105|) we see that when 
B x is very small we have to evaluate B 2+e / 2 (P 2 ; m 2 , M 2 ) , which is of the general form of Eq. l|38l) . With 
a = 2 + e/2 in Eo. H38|l and discarding the terms with positive powers of A(— P 2 , m 2 , M 2 ) we have, 



D 



2+e/2 



D 



1 3 + e 



2' 2 



p) 



3/2+e/2 _ 



3 + e 



[a- 3 - e + (l-a)- 3 - £ ] +0(p- 1 ) 



,-2-e/2 



(114) 



where we have introduced P 2 = — s , p = —(4s 2 ) / A and a = (s + m 2 — M 2 )/ (2 s). It suffices to expand the 
result of Eq. l|114fl in powers of e and to retain the first two terms in the expansion: 



5 2+ e/2 = [' (0) + ^ + O + O (P^)] 



1-a 



j(0) -?V/2_ I J_ _ * 3/2 ( 1 V 31UX2S + 2 



2 r 3 ^— ' A 3 4 
The case b) is quite similar; consider the integral 



A 3 



(115) 
(116) 



Ya = dx x [x{x) — iS\ " = I dx (x — a) [x(x) — iS] a + aB a , 



Since \{ x ) = s i( x ~ a ) 2 + P X ]j the integral in Ea. (|117|) gives: 

^2(1 — a) _|_ ^2(1 — a) 
= ~ 7Z ; V ClB n 



2s(l-a) 



Setting a = 2 + e/2 and expanding in e we get: 



(117) 



(118) 



Y, 



2+e/2 



111 1 e 

2~s I M 2 _ to 2 _ 2 



lnM 2 + l lnm 2 + l 



Af 2 



+ aB 2+e/2 , 



(119) 
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which conclude our discussion. Collecting all pieces together, we get: 



(C 


) r(n 


..,(1 








2B X 


4B x 



i = a.b 



(120) 



.1 

( / dzln K ; m - l) - ft? - + JiiM- (121) 



5.3.2 Evaluation of V C M 

Configuration c) of Fig. is a special case of b). For this case the polynomial V M (see Ea. dlOlf) ') takes the 
form 



V C :m(x) = x 2 M giving a c - M {x) = M . 

1 — x 



(122) 



V c - M (x) has a double zero for x = and the corresponding BST factor is therefore zero. As a consequence 
the diagram acquires an extra infrared pole, confirming the presence of a leading Landau singularity. To 
evaluate this special case, we reconsider Eq. (|93fl . which in the present configuration reads: 



^)r(2 + e)^l|d5 2 ({z})^V 1 [x(l-x)]^/ 2 y e / 2 - 1 (l-z 1 -y)(^^ + Z c;i/ ) ' (128) 



Z C ; M = f3(l- Zl ,l-z 2 ;P 2 ;m 2 ,M 2 ), (124) 

with j3 defined in Ea. (jl2|> . We perform the following set of transformations: y = (1 — Zi)y' , and Z\ — 
1 — z[ , Z2 — 1 — z' 2 , followed by z[ = z'{ z' 2 , obtaining: 



V" = - 



r(2 + e) / dC 4 (x,y,{z}) 



zi y 



z 2 x(l — x) 



e/2 



Z\ 



xyz\ M 2 
l-x 



Z2X(Zl) 



-2-e 



(125) 



where x( z i) = X { z i I P 2 i m ? > Af 2 ) of Eq.©. Integrating over z 2 gives an hypergeometric function and we 
use its properties (see appendix IU|l to get: 

V " = -(t)' 2h r{2 + ,) I "C^y.^x-^-'Hl-xf*"-' z- 1 -'-' 



— / rfC 3 (a;,j/,2) 



r i-- r i + -c x 



-1-26 



(l-a;) 1+e z- e 



xy- 2 -(l- J/ )(Af 2 )- 1 - 3e / 2 X ^^/ 2 (z) - ^-r(2 + £ ) [x(l-,:)]- £ / 2 z 1+e / 2 

x y-^ 2 (1 - „) X - 2 ^{z) 2 F X (2 + e , 1 + § c ; 2 + | c ; - Jgf^ )]. (126) 

In the first term of Ea. (|126p the integrations over x and y are trivial and generate the double pole in e. 
The second term of Ea. (|126[) shows a single pole, hidden in the y integration, which can be extracted as in 
Eq.JBJ. We obtain: 



V£ = T (l - I) r (l + | C ) B(-2 e, 2 + e) - 6, 2) (M 2 )" 1 "^ 2 jf X~ 1+e/2 (2 

r(2 + e) Jdzz 1+ ^ 2 X - 2 - e (z), 



(127) 



2 2 
e 2 + 3e 



K M 2= / dC 3 (x,y,z) 



z (I- x)x(z) +xzM 2 



X 2 (z) (1 - x)x(z) +xyzM 2 



(128) 
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For V C M we perform a Laurent expansion in e which gives: 

1 



V M 
v c,0 



2NP 



r (1 + e) ( 11%' e- 2 + Tl c iX (- 1 + K% M 



(129) 



R«Z = J °, Wzf = ( 1 - ~ InM" ) J" + - Jg - Jl' 



2 



= (| In 2 M 2 - | InM 2 + H C (2) + l) J* + (l - \ InM 2 ) (1 J* - J?) + 1 J 2 - i J, 1 + \ J 2 °, 
where J% is defined by: 



e?z In" z 



In x( z ) 



This class of integrals can be treated by using Ea. l|51[) which, in the present case, reads as follows: 



d z ln fc+1 x(z) 



(130) 



(131) 



After integration by parts, we get: 



7° - 



/i - 



rfc-l _ 
J Q — 



1 



2B x 




1 






L 


1 




2B x 


k ' 



<iz In'' 



Z I 7 



z-X, 



+ ]nz) ln44 + 2(-l)"r(n + l) 



— I dz |4 In™ z lnx(z) + In" 1 z (n 



x(Q) 

z-X 



n > 1, 



In z ] In 



xO) 
x(o) 



[ln X W+lnx(0)]} ; 



dz In*" 1 x(«) pn*(z) + 2fc] - £ x j, 



with L£ given in Ea. l84|) . The same strategy is adopted for the second integral, thus obtaining: 



-=r(i + e ) ( - 2/ ° e - 1 -/ -/ 1 + 2/ 1 



where I* is now defined by: 



1 „ ln fe x(-z) 
dz z In™ z ■ 



X 2 (^) 



This class of integrals can be treated by using Eq. l|52|l which, in the present case, becomes: 
X~ (^) In x(z) = -5- |x~ (2) In H r 2^ ~U dz X ~ ^ ln X ^ 



n> 1, 
(132) 

(133) 
(134) 

(135) 



After integration by parts, we repeat the whole procedure for I*; the BST relation to be used in this second 
step is: 



X^(z) ln fc X (z) = -L ln fe - * ~ ^ S 2 ln fe+1 X («) 



2(fc + l) 



(136) 



The second integration by parts gives: 



~ 452 



1 

dz 



X 



k + 1 



±- ln fc+1 X (z)-(4z-5X x ) \n k X (z)-(Sz-5X x )k\C k - 1 (z) 



X, 



X 



x 



k + 1 



ln K+1 x(l)-(2-Xjfc!£ fe (l) 



X y £ k+1 (Q) + 2 



fc! X x l + £ fc (l) 
2 5 X X(l) 
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where £k{z) — Y)i—i m ' x( z )/'- Finally we evaluate V™ 2 . First the integration over y is performed: 



V* 



o Jo 



dx I dz 



1 - x 



x M 2 x{ z ) X 2 ( z ) 



ml 



xzM 1 



(1 - x)x(z) 



(138) 



In the second term of Eq. i|138|l we transform the argument of the logarithm to get 



Vch = das dz 
o Jo 



l-x 
■M 2 x(z) 



In 1- 



xzM 2 



(l-x) X (z)J x 2 (z) 



In 1 



x z 



M 2 



-IxlM 2 (139) 



where 1% is given by Eg. 11371 1 . The first two terms in Ea. ll39|l are treated according to Eci. (|5T|) and Eci. (|o7]|) 
which in, the present case, become 



In 1 



In 1- 



xO) 
X(z) 



1 
1 



z 1 + 



z 

z-X 



In 1 



x(*0 



x(*0 



In 1- 



z-X, 
2 



<9 Z Li 2 



a z 



X« 
z-X 



2 a 



■ & In 1 - 



x(*0 



where a = -iM 2 /(l — x). After integration by parts, we finally obtain: 



V* 



1 



2B 



x Jo 
X 



dx / dz 



1 - x 



o 



X 



In 1 



iiM 2 

(l-ar)xW 



In 1 



xzM 2 



(1 - x)x{z) 



Li 2 



(i-^)xW 



:M 2 



X 



Collecting all pieces together, we get: 



y<r\=- ( v) r ^ 



M 2 R 



.C(3)+lnM 2 / + / 1 °-/ 1 



2 v -2 



vi 



2 M 2 1 



2M 2 °' 



- vx 



TZ C ' M 
~2M 2 



yC-M 

v 



J° + 2I l + V c M 2 . 



5.4 The V G diagram 

This family of diagrams, corresponding to Fig. [8] can be written as 

1 



7T 4 V/ G = M 2£ 



d n qi d n q 2 



[1] G [2] G [3] G [4] G [5] C 



where we have introduced the following notation: 

[l}o = q 2 + m 2 , [2] B = (qi+ Pi) 2 + m 2 , 



[3] G = (5i - qif +ml, 



[4] = (g 2 +Pi) 2 +' 



[5] G = (g 2 + P) 2 + 



The leading Landau singularity corresponds to the following four solutions: 
1 



P z 



A, 



4m 2 m| 



(m 2 - m 3 + to 4 ) (p l +m 1 - m 2 ) (-p 2 + m 4 - m 5 ) - 4m 2 m 4 (p 1 + p 2 ) 



- V (p 2 - to 2 + to 2 ) (A 12 )!/ 2 + a (p 2 + to 2 - to 2 ) (A 13 ) 1/2 + a V (to 2 - m 2 + to 2 ) (\ 23 )^ 2 
XiXj, Ai = A(to 2 ,to 2 ,to 2 ), A 2 = A(-p 2 , to 2 , to 2 ), A 3 = A(-p 2 , to 2 , to 2 ). 



(140) 



(141) 

(142) 
(143) 



(144) 



(145) 
(146) 

(147) 



where 7717 = ±1. From a study of these conditions it follows that V E is free from QED - like infrared 
divergent configurations. 
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Figure 8: The irreducible two-loop vertex diagrams V G . External momenta are flowing inwards. 

5.5 The V K diagram 

This family of diagrams, depicted in Fig. can be written as 

*' v " = " 2< / [ii.w.w.'m.w.w. - (148) 

where we have introduced the following notation for propagators: 

[l]*=«i+m?, [2} K = ( qi + P) 2 + m 2 2 , [3} K = ( Ql - q 2 ) 2 + ml, 

[A] K = ql+ml [5] K = (q 2 +pi) 2 + ml [6} K = (q 2 + P) 2 + m 2 , (149) 
The leading Landau singularity for V K is given in Eq. (321) of III. Our strategy will be to search for 




Figure 9: The irreducible two- loop vertex diagrams V . External momenta are flowing inwards. 

solutions of the Landau equations where two external momenta are on some mass-shell and the third one 
is unconstrained. Suppose that P 2 is unconstrained, then we must have TO3 = and, moreover, the first 
condition of Eq. (321) of III requires m^ — mi and mg = m%. The second condition will be written as 

a 2 2 2 ( 2 2 2\ I 2 2 2\ I 2 , 2 , 2\ ,1 2 2 / 2 , 2\ 

4 m 3 m 4 p 2 = (m 1 — m 3 — m 4 ) (m 2 — m 3 — m e ) (p 1 + m 4 + m 5 ) — 4 m 3 to 4 (m 5 + m e ) 

+ rji {ml - m\ - ml) A^ 2 - r) 2 {ml -m\- ml) A^ 2 + 773 (p 2 + m 4 + m 2 ) A^ 2 , (150) 

where r\i = ±l,z = 1, . . .3. With m^ = mi and to 6 = m 2 we obtain \i — m\ {m\ — 4 m 2 ) and A 2 = 
m 3 ( m 3 — 4 m 2 ,), so that m 2 . can be factorized on both sides of Ea. l|15U|l . Next, by putting m 3 = we obtain 

m\ (p 2 + ml + ml) + s 4 mi m 2 {pi + m 2 + to 2 ) = 0, (151) 
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from which we recognize the familiar infrared configuration with p\ = —m\, p\ = —m\, and = mi, 
m,Q = 1712,1123 — TO5 = 0. which corresponds to configuration b) in Fig. 1101 If p\ is unconstrained (the 
case with p\ leads to the same configuration because of the symmetry of the diagram), we must have 
m4 = 0. Searching for configurations with two photons, there is only one possibility compatible with the 
standard model, m,2 — 0. Inserting these conditions in the second relation of Eq. (321) of III we obtain 
mi = 1713 = rriQ and p\ — — m§. Then the first relation of Eq. (321) requires P 2 = —m\, leading to 
configuration d) of Fig. In the following we classify the sub-leading singularities. For the sub-leading 
singularity corresponding to ai — 012 — «3 = the reduced diagram is a one-loop three-point function 
where the classification of infrared singular configurations is already known: we obtain three solutions which 
correspond to the two configurations shown in Fig. ^] The first is 7714 = 0, P 2 = — mf, p\ = — m 2 , the 
second is m 5 = 0, p\ — —m 2 , p\ = —rn\, with a third one given by to 6 = 0, P 2 — —m 2 , p\ = —m 2 . 





Pi Pi 



Figure 10: The V K infrared configurations. The photon line represents a general massless particle while dashed and 
solid lines refer to different massive particles. The mass of the three particles in the internal triangle are mi , m,2 and 
m 3 . 



5.6 Evaluation of the V K cases 

In order to evaluate infrared configurations for this family of diagrams we use the parametrization of 
Eq. l|149|l . At first we combine propagators [1] K — [3] K with Feynman parameters xi,X2, 

- 4 = " 2l r (3) / d "^> I «■«■» ra* wW^tw (152) 

where ^-dependent quantities are 

R x = {l-Xi)P- x 2 q 2 , Ql = xi (ml - ml) + x 2 (g| + m\ - m\) + m\ + (1 - x x ) P 2 . (153) 



2G 



Integration over q\ gives 

, e r(l + e/2) 



7T 2 V K = I fi- 



e/2 



dS 2 ({x}) x 2 (l-x 2 ) 



-l-e/2 



d n q 2 



(q 2 + 2P x -q 2 + Ml) l+</2 [4] K [5]* [6] K 



,(154) 



with new auxiliary quantities denned by 
1 — X\ 



Pr 



P = XP, 



Mi 



x(xx ; P 2 ; m% , to?) + x 2 (m§ - to?) 



1 — ara ' £2 (1 - x 2 ) 

Secondly, we combine the remaining propagators with Feynman parameters j/i, i = 1, . . . , 3: it follows 



(155) 



7T 2 v K = i ^ n-^ r( 4+ |) / ds 2 ({x}) [x 2 (i - z 2 )r W2 / ds 3 ({j/}) yf 



X / d n q 2 



y 3 [123] K + (y 2 - y 3 ) [4] K + (1 - j/i) [5]* + (y x - y 2 ) [6] K 



-4-e/2 



(156) 



where [123] K = qf + 2P X • ga + Mj. After performing the ^-integration we obtain the following result: 



V* = - T(2 + e) J dS 2 ({x})[x 2 (l-x 3 )]- 1 -'* J dS 3 ({y}) 



(157) 



where the quadratic form U K is given by 

U K = -[p 2 y 1 -P{y 2 -Xy 3 )+p 1 } 2 + (P 2 -p 2 +m 2 - m 2 )y 1 -(P 2 +m 2 -m 2 )y 2 + {M 2 ~m 2 )y^ (158) 

To proceed in our derivation we introduce X = 1 — X (with X given in Eg. (1155(1 1 and perform a change of 
variables, y 2 = y' 2 + X y 3 , obtaining 



V 



T(2 + e) / dS 2 {{x})[x 2 {l-x 2 )\ 



A rXyi,yi 

''~ e/2 / ^ /_ dy 2 y lm2 , K 

Jo Jo.Xyi 



y t ; K = I dy 3 y\ j2 (A K y 3 + B K ) 
Jo 



-2-e 



0-1 = =■, 

X 



a 2 



yi - V2 
X 



(159) 



where A K and B K are defined by 

A K = Ml -Xmj- (P 2 + ml) X 
B K = -(pa yi-Py 2 +Pi f + (P 2 -p\ + m\- m|) Vl - (P 2 + to 2 - to 2 ) y 2 + p\ + to 2 (160) 

The y^ K functions of Ea. (|159p can be expressed in terms of an hypergeometric function, 

r 2 (l + e/2),. 2 ,_ - e _ e B„ 



y l - K 



(A K B K )-^/ 2 -L- (a, A 2 K )-^I 2 2Fl( 2 + e,l + e -;2 + e -; 

I + e 2 2 a,i A K 



r(2 + e) 

The first term in Ea. (|161|) does not depend on af, therefore, we obtain, 

V k = _ (i\ f d s 2 ({ x }) L a (i _ x ^ 



(161) 



V 2 1 



f T(2 + f) P 1 r x vi-vi -, 

/ dS 2 ({y}) (A K B K )-^I 2 -2 \2 ' d Vl _ dy 2 J 1(B2:K , (162) 
J z + £ Jo Jo.x yi J 



where Ji- K is defined by 



J %K = (a, Aiy 1 ^ 12 2Fl (2 + e , 1 + t ■ 2 + | ; — ^- 



(163) 



Having derived the general result of Eq. l|162f) we proceed with separate evaluations of the four cases of 
Fig.Hni 
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5.6.1 Evaluation ofV* 

For the configuration a) of Fig. QJJwe have: 

V a - K (xi,X 2 ) 



x 2 (1 - x 2 ) ' 



B a ; K = (V2 , Vi ; P 2 , M 2 , m 2 



(164) 



V a . K = -P 2 x { + (P A - m 2 + M 2 ) x x x 2 + m 2 x 2 + (P 2 + m(- mi,) x x - (P 2 + M - m\ + m{) x 2 + m 2 2 . 

After a transformation of variable, y 2 ~ * y\y 2 , in the first integral of Ea. (|162p the j/i-integration can be 
carried out giving raise to the infrared pole. After that we put e = in I\. K and I 2 - K , thus obtaining: 



V K = 



T 2 (l 



dS 2 {{x}) 



dy 2 



V a ;KX(V2) 



-l-e/2 f 1 fXviiVi 

+ V~£ d Vl _ dy 2 C^ 2 

JO JQ,Xyi 



, (165) 



where C%, + B a -K and \ is defined in Eq.©. Henceforth we change variables in the last two 

integrals, y 2 = y\ y' 2 \ the y\ -integration can be performed, giving 



C = ( - ) T 2 ( 1 + - ) / dS 2 ({x}) i-l dy 2 V a -K xM 



7T 



-l-e/2 



1 



In 



1 



x 2 (1 - x 2 ) h ai(y 2 ) 



U2 V a - K 



jr[ Jo ai(y 2 )V a . K 
where we have introduced 

a iW = xW, OL 2 {y 2 ) =x{Vi), bx=X, 
In the last two integrals of Eq. l|166[l we change, once again, variables: 



X. 



x\ = x[ + x 2 , x 2 = 1 — x 2 , x\ x 2 , for i = 1, 
x\ = 1 — x'x, x 2 = 1 — x' 2 , x\ <-> x 2 , for i = 2, 



thus obtaining 



V* 



2\ e 



r 2 i 



<^a;ic + Ja.l-K + Ja,2-K, 



(166) 
(167) 

(168) 
(169) 



where the various ingredients are defined by 



dS 2 ({x})v-*- e / 2 ( Xl ,x 2 ), y. 

i 



r x 2 /xi 

J a ,i-,K = I dS 2 ({x}) I dy 2 



Oii{y 2 ) V a ,i; f 



hi 



1 



dy 2X - 1 - e/2 (V2) 
x 2 (1 - x 1 )a i {y 2 ) 



The quadratic forms introduced in Ea. (|170| l are: 

V a ,i- K {x\,x 2 ) = x l H a ,i x + 2 Kl , t x + L a ,i, 



2/2 V a ,i-K 



* = 1,2; 



(170) 



(171) 



V a ,i; K = M 2 xl + {P 2 - M 2 + m 2 ) xi x 2 - P 2 xl + (ml -m\- M 2 ) X\ + (m\ -mj + M 2 - m 2 ) x 2 + m 2 , 
V a ,2-,K = m 2 x\ + (P 2 - m 2 + M 2 ) X\ x 2 - P 2 x 2 2 + (m\ -m\- m 2 ) X\ + (ml - m\ + m 2 - M 2 ) x 2 + m\. 

To compute X a - K and y K , defined in Eq. (|169|l , we use a BST relation. Given V a - K = x t H a x + 2 K l a x + L a 
we introduce X®. K — 1, X\. K = and also 

B a , K =L a -K t a B- 1 K a , X^ K = -{K l a H-^ i = 1, 2; (172) 
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we obtain for X a;K = X\. K + X 2 . K e 
(_l)»+i 



■yn 

a ' K ~ 2 n n\B a . 



2 f dS 2 In"" 1 V a , K ( In V a:K +n) (X« :K - X^) f dx In" V£. K 

J V 7 *=o 70 



where the V*- K quadratic forms are 

K;k = V a , K (l,x), V?. K = V U . K (X,X), K;k = V a;K ( X) 0). 

Similarly, for y K , we get 



y K = yi + y\ e, 



y n K 



(~i)" +i 

2" n! R, 



<2y In" ^(y) lnx(y) + 2n 



(173) 



(174) 



(175) 



where BST factor and co-factor are defined in Ea. (jl()|l and in Eq.JHU). To compute J a ,i-K we use Ea. (157(1 
which in the present case reads: 



I ]n ( 1 _ A **K 



i 



i hl ( 1 _ A <>*« 



V a . 



3=1 



5 fe)}- 



(176) 



where ^4 a .i:K = — x 2 (1 — #i) ct{/y 2 and BST factor and co-factors are 

S a; /f = L 0l i-if* )1 i?'~ 1 Jf 0)1 = La^-Kl^H'^Ka^ = La-KlHa 1 K a , 
After integration by parts we obtain 



1 



B a -K Jo 



dy 2 



f xi 1 
ax\ I dx 2 — 7 — r 

Jx lV2 a i{VV 



1 — Xi 



•X'2 



X M;iC = -if*,lf Q - 1 . (177) 



ln(l + rja.i-K) + Li 2 (-%,•*;*:) 



dx 



a.i:K aA:K T . / l \ y2 ^a.z'l-fC a.z;K T . / 2 \ 
Ll2 (-^o,i; K j —7-I Ll 2 (~»?o,*;icj 



a»(y2) 



where the 77-functions are defined by 

_ X 2 (1 - Xi) aj(y 2 ) 
?/2 Va,i;K(^l,^2j 

Collecting all pieces together, we get: 



x(l-x) ai{y 2 ) 
y 2 V ati;K (x,x) 



}■ 



2 _x(l-x) ati(y 2 ) 



Va }i . K (x,xy 2 ) 



r(l + e) (i^T + C) 



■^K^ _ X a; K yl ~ Ja,l;K ~ J a ,2;K- 



(178) 
(179) 

(180) 
(181) 



In the following section we will consider the b) configuration of Fig. ^3 



5.6.2 Evaluation of V h * 



Configuration b) of Fig. ^| is a special case of a). Here the polynomial V takes the form ((3 is defined in 
Eq.©) 

V b;K {xi,x 2 ) = (1 - x x , 1 - x 2 ; P 2 ; m 2 , M 2 ) (182) 
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which has a double zero at x\ — x 2 = 1. The BST factor for the polynomial Vb-. K is zero, revealing again 
the presence of a new infrared pole. For this reason we cannot put e — from the beginning; instead we 
reexamine Eq. l|159| ) which, in the present configuration, reads: 



V, 



, A ' = -(^) £ r(2 + e) J dS 2 ({x})[x 2 (l-x 2 )}- 1 ^ 2 J\ yi 

/3(1 -xi,l- x 2 ) 



(0,X, V1 ) 



dy 2 dy 3 y 3 



x 2 (1 - x 2 ) 



■2/3 + /%2,2/l) 



In the first (second) integral of Ea. l|183|l we change variables according to: 

xi^xi+x 2 , x 2 ^l-x 2 , xi^xix 2 (xi— >1 — Xi, x 2 ^l-x 2 , x\ 
thus obtaining 

[x 2 (1 - x 2 )]- e/2 ^ f h * f d * e/ 2 [X 2 ai(xi) 
■ y. dy 2 dy 3 y 3 ' 

,; =1 J a, JO 



Xi X 2 



V h=-[^) T(2 + e)ldCs({x}, yi ) 



1 — X 2 



2/2 



1 — Xi 

2/i - 2/2 



2/3 + ^(2/2,2/i) 



ai=0, h = xiy!, dx = — , a 2 = (1 - £1)2/1, 62=2/1, da = 



(183) 



(184) 



-2-6 



(185) 



where the functions a.i are defined in Ea. (|167p . We also transform the integration domain of y 2 and y 3 in 
the first(second) integral of Ea. fl83|> according to: 



2/3 -» 2/2 2/3 M, 2/2-^2/12/2, (j/3 (2/1 - 3/2) 2/3/^1, 2/2 -» 2/1 (1 - 2/2)) ■ 



(186) 



After this transformation we obtain 



VJ, K = y T(2 + e) ^ I dC 3 ( Xl ,x 2l y 3 ) f 



^2/2 



sci (1 - x 2 ) y 2 



l+c/2 



X 2 



yb,i- K = I dyiy 1 t/2 x 2 y 2 y 3 ai{xi) + xi(l ~ x 2 )yicti(y 2 ) 
Jo 1 

yb,i\K is then expressed in terms of an hypergeometric function: 

y b , t ; K = [*a 2/2 2/3 Oi(a:i) ]~ 2 ~ e 2^1(2 + e , 1 - | ; 2 - | 

= S(l - 1 1 + I e) [x 2 y 2 y 3 a t ( Xl ) J- 1 ^/ 2 [a* (1 - x 2 ) ^(2/2) T 1+£/2 

- w^-[xi(l-x 2 )a i ( y2 )]- 2 -' 2 F 1 (2 + e, l + |e; 2 + |e; , 
2 + 3e 2 2 xi (1 - x 2 ) a* (3/2) 

Substituting this partial result in the expression for V b K and setting e = where possible, we get: 

where we have defined 
2 ,, 

e ,,-1-e -l-3/2e. 



-2-e 



(187) 



x\ (1 - £2)0^(2/2) ^ 
x 2 2/2 2/3 aifal) 

+e/2 

x 2 2/2 2/3 Oi(a;i) 



(188) 



(189) 



dy 2 x l x 2 (1-^2) 2/2 2/3 (^l)", (2/2) 



y M=E / ^3(^2,2/3) / 

^50 = / ^3(^1,^2,2/3) / ^2/2 ^ 2 ; [^22/2 2/3 + xi (1 - z 2 ) 0^(2/2) T 1 - 

Jo «i (2/2) 



(190) 
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After some straightforward integration the function of Eq. I|19U|) become 



v,; 



1 r(l-2e)r(l + e) 



2e 2 



r(l-e) 



2 pi px 

V dx dy x e y~ e a; 1 " 3 ^" {x) a; 1+e/2 (y) 
._, Jo Jo 



2 p p Xl 

4=1 7 Jo 



ln(l + T] b}i . K ) 
ai(xi)ai(y) x 2 " 



VbS;K 



x 2 yaj(x 1 ) 
X\ (1 - x 2 )ai{y) 



After a Laurent expansion in e, V b K has the following form 



Ke = l E ^ e "' 



n=-2 



K 



' fdx X -\x)\ 2 = [j t ^ = - j ° + y: to - •»)' 

9 

= - 1 US? + \ Jl j°o + 2 cm [j °] 2 + i E [- 2 J n ( j ) 

i=l 

+ (<) + - 3 ^o°(i) + ^o 1 (0 + 3 Jo'iW - (<) 
In Ea. i|193[l we have introduced and Jnm(i) defined by: 



J* = \\ z ln « , !^Cf) Jfe , (i) = /"^ A ln „ x ln ™ In* «,(*) In' 1 
7 Jo h 



a>i(x) ctiiy) 



(191) 



(192) 



(193) 



(194) 



J integrals of the first kind have already been considered in section 15.3.21 (see Ea. (|132|) '). For J integrals of 
the second kind it is easy to verify that the following property holds: 



fe 2 





(195) 



The general analysis of J is reported in appendix [DJ here it is enough to say that for J we are always able 
to find a smooth integral representation. Finally, note that the residue of the double pole is exactly one half 
of the square of the residue for the single pole in the one- loop case. For V b K we apply Ea. (|57|l : 



1 Infl- by 



1 



y- Zi 



d y Li 2 



by 



a i{y) ^ OLi(y)J B x ]\ 2y I V ai(y)J 2 \^i(y) , 

where b = —x 2 Oi(xi)/xi/(l — x 2 ), Z\ — X x and Z 2 = X x After integrating by parts, we obtain: 



V b,0 



i 2 r c Xi i \t z \ 

— r Y, dC 3 {{x}) / dy — 1 + — ln(l + %4;A -)+Li 2 (-7 7M;fc ) 

2B xl~{J Jo x 2 a l (xi) [\ y) 



Collecting all pieces together, we get: 



7) r (! + e ) ( ^ V -2 + \ V -i + V o ' K 



v b j 2 K = \n% 



V-l "2 K -l 



(196) 



(197) 

(198) 
(199) 
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5.6.3 Evaluation of V C K 
For configuration c) of Fig. ^]we have 
A _ V c . K (x 1 ,x 2 ) 



X 2 (1 - X 2 ) ' 

B clK {yi,y 2 ) = P (2/1 - 2/2 , 1 - V2 ; P 2 ; to 2 , M 2 ) 



V r c;K (xi,a;2) = 5 -M 2 ; mf , to 2 ) + x 2 (to 2 - to 2 ) 



(200) 



For this configuration we transform variables according to y\ = y[ + y 2 and 7/2 = 1 — y'2 m the first integral 
of Eq. Hl(j2f> which will have the same form as the similar integral in case a) and we obtain: 

m 2 V 1 



v c = [yj - r a (1 + -J * C;K y K + j cA . k + j ca , K , (201) 

where the various components are defined by 

X c ,k = J dS 2 {{x})V-i- e/2 ( Xl ,x 2 ), y« = J Q dy2X~ 1 - e/2 (y2), 

//■l rX yi,yi 

dS 2 ({x}) / dyi / _ dj/aC^ajciK) ^c;* = «i ^ c; k + B C . K . 
JO J0,Xyi 

and x( x ) is defined in Eq.©. y K has already been computed in Ea. l|175[) . while for computing X c;K we can 
use the following BST relations: 



(202) 



(203) 



1 



TO3 — m{ [i+l 



8 V^ +1 

U X2 v c;K ' 



r CI-ftT 



/i + 1 



A(M 2 ,to 2 ,to|) 



X C ;k — 



M 2 — TO 2 + TO 2 



4M 2 ' 2M 2 
After intergration by parts, we obtain the following two results: 



(xi - X C . K ) d Xl + x 2 d X2 



X C ;K — 1 ~ X C:I 



yfi+l 
v c;K ! 

(204) 
(205) 



1 



4* -T^J if ™l^™3, * c;/c 



2 6, 



^ci-T^i if b o;K^ 0- (206) 



Afi" - / dx 



ln n V c;K (x,x)-\n n V c;K {x,0) 



(207) 



A- 2 " = / dS 2 In 



n - 1 V c; *(3mV c;K + 2n)- / dx\x c , K \n n V c)K {l,x) + {x + X c , K )\n n V c , K {x,x)\. (208) 



Therefore, our result is also valid when mi = to 3 or A(Af 2 , to 2 , to 2 .) = 0, but not when the two conditions 
occur simultaneously. In the standard model this can actually happen only when mi = to 3 = M and TO2 = 
which corresponds to configuration d) to be analyzed in the next subsection. 
In order to compute J Cj i ;JC and J c . 2 - K , we change variables according to: 



xi = x[ +x 2 , x 2 = l- x' 2 , xi <-> x 2 , y 2 = yi y' 2 , for i = 1, 

xi = l-Xi, x 2 = l-x' 2 , Xi ^ x 2 , y 2 = y x (1 - y' 2 ), for i = 2, 



thus obtaining 



/■ pi rx 2 /xi ^ 

Jc,t- K = / dS 2 {{x}) / dt/i / dy 2 — — 

J JO JO 



x 2 (1 - £1)2/1 



yi 2/2 + x 2 (1 - xi) -B c ,i; K ' 



(209) 



(210) 
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where we have introduced new quadratic forms Vc.i-^K- 



V c ,i;k = x(l -xi+x 2 ; -M 2 ; mj , mf) + (1 - x\) (m§ - mf), B c ,i;jc = B C]K (yi, 2/12/2) ( 211 ) 



K,2;i< = x(l - ^2 ; -M 2 ; , mf) + (1 - Xi) (m% - mf), 



^,2:^ = B c , K (yi, yi(l - y 2 )). 



For V c ,i;k and V^,2; a- two BST relations are available. Looking at the expression for J c ,%;k, we immediately 
see that the required BST relations are of the form of Eq. (|56|l . We have 



1 



1 A 1 



Vc,i\K A 



Vi,iA + V l , 1 In 1 



.4 



where the explicit form of 7\j is 



1 



Vc,i;K A 



(l+V h2 )A + V h2 In 1- 



A 



(212) 



Pi, 



-d Xl - d X2 , Pi. 2 = (1 - xi) d Xl + - (1 -xx - x 2 + X C , K ) d X2 , 



V 2 a = (1 - xO + - (1 - x 2 - X C;JC ) <9 K2 . 



7\i = -d xi , 

In our case it is A = — x 2 (1 — Xx) B c ^ K /(yiy2), and we can define djj such that 



Vi.\ A = a 4 ,i 



2/i 2/2 



ai i B 



2,2 ^cA\K 



2 2/i 2/2 



For dy we obtain: 



ax,i = 1 - xi - X2, ai,2 = (xi + x 2 - X C;A - - 1)(1 - xi) 

02,1 = -%2, 0-2,2 = - (1 - X 2 - X c;K ) (1 - Xi). 



After integration by parts the J c ,i- K can be expressed as 
d 



Z nl 2 £ J tl 



if mi ^ m 3 , 



2b r . 



1 d 



(213) 
(214) 

(215) 



(',2 



(216) 



if b c . K ^0. (217) 



J i,j=- I dS 2 {{x}) Qij / dyi / dy 



X2/X1 2/i 



72 



2/1 2/2 VJ.,i;jr + X2 (1 - Xi) -B c ,i 

X2/Xl 2/i , A , x 2 {l-x 1 )B c . l , 



Ji,i = / ^2({x}) / dift / ' dy 2 /;,.,. In ( 1 

Jt d - j dc 3 (x {y} ) cy ^ m f i + f^^T^) ■ 

J 7 -B C:i;K V 2/i 2/2 K, l;J f(x,x) / 



2/i 2/2 V" c , 



J*. = / dC 3 (x {y}) dkj -p- In f 1 + ^ ''/^ ' '■,). (218) 
,J J B c>i - K V 2/i V" c , i;K (x,xy2) y 

where a^j is defined in Ea. l216|l . while the b ■ ■ ■ d coefficients are: 



(219) 



1,1 


= o, 


Cl,l 


= o, 


di,i 


= -(1-2/2), 


2,1 


= 0, 


c 2 ,i 


= -1, 


<h,i 


= 2/2, 


1,2 


= -3, 


Cl,2 


= i — X c - K , 


di,2 


= l- x-2y 2 + xy 2 + X C . K , 


2.2 


= -3, 


C2,2 


= 1 - x + X C . K , 


^2,2 


= l- 2y 2 +xy 2 - X C . K . 
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All the quantities in Ea. i|218| ) but J?- are stable under numerical integration. For this integral we perform 
back the transformations of Eq. l|2U9[l and get: 



J? d = J dS 2 {{x})a' id J? 



where the new coefficients a' „• are 



a' 11 = Xx-2x 2 , a' 12 = (2x 2 - x\ + X C , K ) x 2 , a' 21 = x l -X2, 
The new integrals Jf are defined by 



Ji = d Vl dy 2 
Jo Jo 



1 



J 2 ° = dyx dy : 



Vi 



2/2 V c -K + x 2 (xi - x 2 ) B c 

1 



ix yi (yi-y2)V c . K + x 2 (i-x 1 )B c . K ' 

We transform according to y\ = y'i + 2/2 > 2/2 = 1 — 2/2 , 2/1 *-> 2/2 , obtaining 



(220) 

a 2,2 = ( X c:j< -Xl +X 2 )X 2 . (221) 

(222) 
(223) 



Ji = dyi _ dy 2 -77- 

Jx J(l-ji)X/X \ L 



1 



2/i) V c -K + x 2 (xi - x 2 ) /3(y 2 , yi) 



ja _ 
J 2 — 



JO 



Hi 



dyx I dy 2 - dyi 



x 



yi 



dy 2 



2/2 V C . K + X2 (1 - xx) (3(y 2 , yi) ' 



(224) 



(225) 



l(X- yi )X/X 

where the quadratic form (3 given in Eq. l|12|l . For J 2 we further transform y 2 — yiy' 2 ; as a consequence the 
2/i integration becomes trivial giving: 



J 2 Q = / dy 



1 



o X2(i-xx)x(y) 



hii 1 



x 2 (i - xi) 2 x(y) 



y[l-xx + {xx - x 2 ) y] V c;l 



(226) 



The x function is defined in Eq.@ . For J" we can use the fact that the B factor of the polynomial (3 is 
zero. In this way we can use Ea. (|58|l : 



V - A 



Vo+V t 1 d x )A = ~V t 1 d x ln(l- 



(227) 



To derive Jf we replace 

V^p(y 2 ,yx), 
and obtain the equation 



A 



(i-yi)Vc; K 

X 2 {Xx ~ X 2 ) ' 



Vo+V[d x )A = -{2-yx) 



'P t xd x ^-(yxd yi +y 2 d y2 ), (228) 



Vc-l 



x 2 (xx - x 2 ) ' 

Substituting in the expression for Jf , Ea. (|224(l . after integration by parts, we get the final result for J", 



(229) 



Jf = 



1 



V C , K 
1 



dyx 



x 



Vi 



dy 2 — Its <! I 



(X-yi)X/X ( 2 -2/i) 



(i-yi)v c , K 



x 2 (xx - x 2 )f3(y 2 ,yx) 



X 



dy =- In i 1 



V c - K Jo l + Xy 
Collecting all pieces together, we get: 



yxv c ., t 



x 2 (xx - x 2 ) (3(X y, 1 - X y) 



{vF]=-(^)^^ + ^{-V-x K + Vt 



(230) 



(231) 



34 



yC-,K 



1 



2 d 



o 3 - /, n 

1 



— mt 



i=\ A=a 
2 d 



2 6 c; k 



x 21 y 1 



yC-K _ 



2b r _ 



5. 6. 4 Evaluation of Vf 

Configuration d) of Fig. ^| is a special case of c). In fact the polynomial V K takes the form V d . K — M 2 x\, 
which has a double zero for x\ = 0. The BST factor for the polynomial V d - tK is also zero revealing the 
presence of a new infrared pole; we are not allowed to put e = but reexamine instead Eq. l|159[l which now 
is: 



v d =~ r(2 + e) / dS 2 ({x})[x 2 (l~x 2 )} 1+ ^ 2 / dm _ dy 2 yw2;d;K, 

JO JQ,Xyi 

-2-e 



yd,i- K = / dy 3 y 3 



, e / 2 



(I 



M x(y 3 + x 2 (l~x 2 )B d;[ 



2 + e 



[X2{l-x 2 )]- 2 ^ ^ 



.l+«/2 
2^1 



d:K 



aiM 2 x\ 



x 2 (l - x 2 )B d . K 



(233) 



where a x = y 2 /X, a 2 = (yi - y 2 )/X and B d . K (yi,y 2 ) = f}(y x - y 2 , 1 - y 2 ; P 2 ; m 2 , M 2 ), with j3 defined in 
Eq.lO. We have also defined the shorthand notation 2 F\ (x) = 2 F\{2 + e , 1 + e/2 ; 2 + e/2 ; x). Now we 
split the integration region, [Xy\ , y{\ = [0, z/i] G [0, Xy{\, and use the properties of the hypergeometric 
functions to obtain 



(234) 



Jd,0;K 
Jd,l;K 



T(2 + e) 
l + e/2 



dS 2 ({a;})/d52(M) 



2/i - 2/2 



-| l+e/2 
2 F\ 



M 2 x{ {yi - y 2 ) 
x 2 (l - xi)B d . K 



1 

M 2 



^({x}) / ds/i / c?2/ 2 o 







1 







^,2;a- = --^2 J dS 2 ({x}) dy x J dy 2 ^ \, 



M 2 a; 2 2/2 + 2:2 (a;i - £2) -Bd ; *r ' 
1 



^2 (1 - Xi) 

2 M 2 xl( yi -y 2 )+x 2 {l- Xl )B d ., 



(235) 



We make the transformations x 2 — x\ x' 2 , y\ — y[ + y 2 , y 2 = 1 — y 2 , , yx <-> 2/2, obtaining 



Jd,0;K — — 
Jd,Ui 



| dCaCM) | d5 2 ({2/})xr /2 



2/2 



£2(1 - xi)(3 2 (y 2 ,yi) 
x 2 (1 - X2) 



l+e/2 

2^1 I 



M 2 a;i y 2 



£2(1 - x 1 )P(y 2 ,y 1 )J ' 



^ J d c 2 ({x }) jj yi J {i _ yi)xt _, dy 2 M2{1 _ yi) + X2 (1 _ X2) ^ yi) •. 

1 ryi 



Jd,2: K = -^H [ dC 2 ({x}) [ dy x [ dy 2 ^^ — — — ■ ^ — 

M J Jx> J(x- yi )X'/x' xx M 2 xxy2 + x 2 {l-xx)f3{y 2 ,yi) 

where is defined in Ea. (|12|l and 



, (236) 



X' 



1 — Xx 
1 — Xx x 2 



x' = i-x'- Xl{1 - X2) 



1 - Xx X 2 



(237) 
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For Jd,i-K (Jd,2;ie) we a Pply the same technique already used in the previous subsection to treat Jf (J 2 ) 
We obtain 



Jd,i- K = -^j J dx 2 J dS 2 {y} 



^ / dC 3 ({x},y) 



(2-yi) 2 

x 2 (1 - x 2 ) X' 



(l-yi)M 2 



l + X y 



hW 1 



1 — 2/i ) "~ I ^2 (1 - x 2 )/3(y 2 ,yi) 
x\ y M 2 | 



Jd.2:K 



M 2 



dC 3 ({x},y) 



Xl 



x(y) 



hit l 



^(l-^^/JCX'^l-X y) J 

x 2 (1 - xi) (1 - x 2 )yx(y) 

[1 - xi + xi (1 - x 2 ) y] [M 2 xi y + x 2 (1 - ari) x{v)] 



(238) 



To extract the infrared poles from Jd,o.K of Ea. (|236|l we first transform according to y 2 = yiy' 2 and then 
perform a sector decomposition of the unit square in xi and yy. [0, l] 2 = [0,1] ® [0,xi] © [0,1] ® [0, yi]. 
Each sector is then mapped back into the unit square giving Jd,o.K = Je,i + Je,2, with 



j eA = -(M 2 )- 1 - £ / 2 r 2 (i + e /2) dc 3 (x,{y}) x 1 1 - 2e y 1 1 - e x(y2)- 1 - e/2 



(239) 



J, 



e>2 



T(2 + e) 
l + e/2 

r(2 + £ ) 

l + c/2 



dCidxj^y^x^y, 1 - 21 



x 2 (l-xi)]^ 2 

2-ri 



2/2 Af 4 



^2 (1 - xi) yi xivz) 



in 



x 2 (l - xiyi)x 2 {y2) 



M 2 y 2 

l+e/2 / , r 2 

2 -Pi 



M 2 xiy 2 



X2{1 ~ xiyi)x(y 2 ) 



For J £i i we have used the properties of the hypergeometric functions of related arguments (see app.EJ- In 
the first term of J £i i, the integration over xi and yi is trivially performed giving the double infrared pole. 
The two terms containing the hypergeometric function are also divergent for e — > 0. The corresponding pole 
can be extracted as in Eq.®. Setting e = when possible, we get: 



Jd,0;t 
Jd,0:t 



1 i_ (m 2 )- 1 -/ 2 r 2 (1 + ~)y K - *g±g I <ic , ( { r ).{,,)) ./,,,,„ . 



2 e 



l + e/2 



1 

Ye 



-e/2 



X 2 

y 2 M* 



l+e/2 
2-Fl 



x 2 yix(y2) 

M 2 y 2 



2/2 



X2X 2 (2/2) 



l+c/2 
2 -Pi 



M 2 xi y 2 
X2 X{V2) 



1 — Xl 



x 2 1 

M 2 xi x 2 (1 - xi) yi x(ya) + M 2 y 2 



y2 1 



xO/2) yi x 2 (1 - xi yi) x(y2) + M 2 xi y 2 



(240) 



y K has already been considered in section 15.6.11 The coefficients of its expansion y K = y^ + y 2 e + y% e 2 
are given in Eci. (|175|) . Now we expand the hypergeometric functions around e = according to Eq. 
and perform everywhere the xi and yi integrations. Using properties of the dilogarithm we obtain: 



- r(l + e) f 1 dxdyc 1 1 



Li 2 (Vd; K ) - lnry rf;Jf ln(l - rj d - K ) + 



1 x 1 
-ln-+lnx(y) lnr)d; K 

2 y J 

xx(y) 



Vd-K 



yM 2 



(241) 



The result of Eq. (|241l) can be cast in a form suited for numerical evaluation, by introducing the integrals J% 
already treated in section 15. 3. 21 (see Eq. l|132|l ). Combining all terms together we finally write Jd,o-,K as 



d\i< 



(242) 



" 2 " 2M 2 ' 



ybjK 



2M 2 
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Vr 



b.K 



i r i 



2 M 2 



\ I (In 2 M 2 - 4 In M 2 - 10 C(2) - 8) J ° + \ (5 In M 2 + 6) Jq 1 - i (In M 2 + 2) J° 
^8 4 2 



2 Jl 8 J " 2 < 



2 t 2 / da; dy — - 

x(y) 



L12 + lnr?d ; /< ln(l - rj d . K ) > - J d ,i- K - Jd,2, K - (243) 



5.7 The l/ fl diagram 

Finally, wc consider the V H -family given in Fig. 1111 which is representables as 

1 



7T 4 V H = M 2e / d n qi d n q 2 



[1]„[2] H [3] H [4] H [5] H [6] H ' 



with propagators 



[l]„ = g 2 + m 2 , [2} H = { qi - p 2 ) 2 + ml [3] a = ( qi - q 2 + Pl ) 2 + ■ 
[4]h = (<?i - 92 -P2) 2 +ml, [5]„ = g| + m 2 , [6]„ = (q 2 - pi) 2 + ' 



i 6 - 



(244) 

(245) 
(246) 




Figure 11: The irreducible two-loop vertex diagrams V H . External momenta are flowing inwards. 



This diagram has only one infrared configuration (shown in Fig. 112(1 which corresponds to mi = m§ = 0, 



and p\ 



2 2 
-m p 2 



i\ = —M 2 At first we combine propagators [1] H — [4] H with 




P2 



Pi 



Figure 12: The V H infrared configurations. The photon line represents a general massless particle while dashed and 
solid lines refer to different massive particles. 
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Feynman parameters Xi,X2,%3, 



^V H =/i' e r(4) / d n qi d n q 2 / d5 3 ({a:}) 



1 



[5] H [6] H (ql + 2R x -q 1 + Qiy> 
where ^-dependent quantities are 

Rx = -X2 Q2 - X1P2+ x 3 P, Ql = x 2 ql + 2 x 2 Vi ■ <Z2 - 2 x 3 P ■ q 2 
Integration over q± gives 



with new auxiliary quantities defined by 
1 



-2-e/2 



fe 2 + 2 P a - g 2 + M|)2+e/2 [ 5 } H [6] H 



Pr 



X2 (1 ~ X2) 



X2 (1 - a^i)P2 - 2:3 (1 - 352) -P 



Mi 



x 2 (1 - X2) 



Xl P2- X 3 P 



(247) 



(248) 



, (249) 



(250) 



Secondly, we combine the remaining propagators with Feynman parameters y\ , 1/2 ■ it follows 
,26 , ^ r r 



^ = *^7a r ( 4 + I) / [X2 (1 " x *T*~ e/2 ^ 2 (M) 



>2/ 2 1+e/2 



(, ?2 2 +2i?,. (Zl +Q2 ) 4 + ,/2- 



(251) 

where y-dependent quantities are R y = y 2 Px ~ (Hi ~ 2/2) Pi and Q 2 = j/ 2 Mj. Integration over q 2 followed 
by a transformation y 2 = y± y' 2 , gives 

V* = - T{2 + e)JdS 3 {{x})JdC2({v})[x2{l-X2)]- 2 -^ 2 vr /2 vl +t/ \A H + B HVl )- 2 -\ (252) 



where A H and B H are 

f-2„. 2/2 



4* = M£y 2 - 



x 2 (1-2:2) 



- r 1 2 -f 

P(x3,xi), B H = -\P x y2-pi{l-y2)\ i = &{ 



X2~X 3 
X2 



2/2, 1- 



2:1— X2 

1 — X2 



2/2 (253) 



and the quadratic /3 is defined in Ea. (jl^|> . To proceed in the evaluation of the diagram we perform the 
following changes of variables 



1 — X2 

x 3 ^x 2 x 3 , y 2 — > —2/2, x 2 ^xix 2 , x\ 



Xl - x 2 



V2^xiy 2 , (254) 



and derive 



1 - (1 - xi)x 2 ' 

V H = -(^j r(2 + e) J dC 4 ({x},y2)x- f - /2 [l-(l-xi)x2Y[x2(l-x 2 )}- 1 - e/2 y 1 2 +e/2 y H 



y f i = j dyiy 1 e/2 [An 'rB„ >n 
where A H , B H factors are 



" 2 ~" = ^; 2 -^(^.i-f;2-I;-£> 



, xiy 2 _ ( < 

•Ah = ^ zX(X2X 3 ), 

X2 (1- X2) 



B H = 0(1 - (1 - 353) 2/2, 1 - si y 2 ). 



(255) 



(256) 



The quadratic form % is defined in Ea. pi|) . Using properties of the hypergeometric function and setting 
e = whenever possible, we obtain: 



V H = - 



v) r ( 1 -|) r ( 1 + l e )/ d ^{x},y,)B-^ 2 V" 



(257) 
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h _ -1-2 e 



V H = x\ 



1 - (1 - xi)x 2 



!J2 



V2 



X\ V2X(X2 X 3 ) + X 2 (1 - X 2 )B H 



(258) 



The first term in Ea. (|258l) is an integral of the type shown in Eq.© and the V" diagram can be rewritten 
as 



r i 



ru + -e 



i 
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J X Jy J 



(259) 



J x = / dx 2 x 2 (l - x 2 ) 2t x~ 1 ~ 3€/2 (x 2 x 3 ), 
Jo 

1 



J = / dC 4 ({x} , y 2 ) 



X\ X(X2 X 3 ) 



Jy= [ dy 2 yrx- 1+t/2 ^'{l-x 3 )y 2 ), (260) 
Jo 

x 2 (1 - x 2 ) 



1 



x\ y 2 x(x 2 x 3 ) + x 2 (1 - x 2 ) B H B° H 



(261) 



Here we define B° H = B H (x\ = 0) = x(l — (1 — £3) 2/2)- The quadratic forms corresponding to J x and J y can 
be treated according to: 



X^{X2 X 3 ) 



1 



^(1 -(l-aia) aa) = _ 



2(A*+1) 
1 

2(^+1) 



^2 



2/2 



X 3 

X 



X 



1 - x 3 



J l)2 



X» +1 (x 2 x 3 ) 

X« +1 (l-(l-x 3 )y 2 ). (262) 



After integration by parts we write J xy as 



Jx = Yb~ l dX2 [I ° + '* e] Jv = Yb~ l! dm [I °y + J y (] 



(263) 



I° x = \nx{x 2 x 3 ) + ^ In ^ - ln^a) + 2 



H = 



-- In x(x2X 3 ) + (\nx 2 + 2 ln(l - x 2 )) lnx(a; 2 x 3 ) 
4 x 2 x 3 



X x ln X( x 2 x s) 



X(0) 



I—X2 



X 



x 
x 3 



In 



X{x 2 x 3 ) 3 ^ x 



X(a=3) 



4 £3 



X(^3) 

x(o) 



_ _ 1 3 2 _ 

lnx(x 3 ) + lnx(0)J + - In x&z) - 6 



I° y = lnx(l - (1 - 13) y) + 7—^ — In 4^ ~ Inx^s) + 2 



1-X3 x(i) 



1 



4 = 7 ln2 x(l - (1 - 2:3)2/) - Iny lnx(l - (1 - 13)?/) + 



4 

1 X 



£3) 



In 



X(l - (1 - £3)2/) 

x(i) 



x 



In 



X(x 3 ) 



4 1 - S3 



lnx(x 3 ) + In x(l)J - j hi 2 x(^ 3 ) + 2 



(264) 



The computation of J, Eas. l259fl l|261|) . is actually more involved. First of all we transform variables 
according to x\ = y 2 X\. For the term proportional to B H we write l/x\ = (1 — xi) 2 jx\ + (2 — x\). Now J 
is given by: 



J= I dC 3 ({x}) [ dy 2 =\ — 



ar 2 (1 - s 2 ) (1 - si) 2 _ J_ 

XlX + %2 (1 -X2)0 P x , 



We have introduced shorthand notations: 

X = x(x 2 x 3 ), (3 = (3(1- {l~x 3 )y 2 , 1-xi), 



a: 2 (1 - x 2 ) (2 - xi) 
xi X + X 2 (1 - x 2 )/3 



=x(l- (1-2:3)2/2). 



(265) 



(266) 
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Both terms of Ea. l|265[l are split, i.e. 



J — Jll + J\2 + J2I + J'. 



f 1 

Jll = / dy 2 — - 



11 T J 12 T J 21 T J22, 

(1 -JQi) 2 _ J_ 



dC 3 ({a:})4, 



J 12 — 



1 1 (1-X!) 2 

at/2 = — — =j 

xx P %1 X + X2 (1 - ^2) /? 



J 21 



1 I-X3 £ 2 (1 -x 2 ) (2 - £1) 
ay 2 — — — — — = , 

xi X xi x + x 2 (1 - x 2 )P 



J'22 = fdy* X 4 ^-^n 2 "?^ ^ 

J xi X xx x + x 2 (1 - x 2 ) p 



In order to compute Jn, we use Ea. (j262|) to increase the power of x by one unit, while for j3 and (3 X we use: 



(3 (3 X 



2 5, 



1)2 



l-{l-x x )X 
I-X3 



^d y2 \n]3- (y 2 - — - 2L) ln^ 



(268) 



The presence of (1 — xi) 2 is crucial in preventing spurious singularities at x x = 1. This is the reason of our 
splitting \/x\ = (1 — xi) 2 jx\ + (2 — x x ). After integration by parts we get: 



Ji 



4B 2 

x 



dC 4 ({x} , y 2 ) I° x hi 



111 



1 — x x 

Xl 



0o 



02 



^3 - X. 



hi 



2(1 



CCl(l-X 3 ) /3 2 /?i l-*3 /? 



(269) 
(270) 



where I x is given in Ea. (|264|) and the new quadratic forms are given by 



Pa = Pi x 3 + (1 - xi) (l-x 3 )y, 1 - xi ), /3 00 = x(x 3 + {1 - xi) (1 - x 3 )y ), 
^=£(x 3 , 1-aJi), ^io=x(^3), (271) 

/3 2 = /3(l-(l-s 3 )a:i, 1-n), /3 20 = x( 1 - (1 - a*) x x ), 

It can be easily seen that the logarithms vanish when the denominator of the corresponding factor is zero. 
The result is then smooth enough to be integrated numerically. 

For J12, we use Ea. l|56|) which in the present case reads as follows: 



1 (1 -jri) 2 
P I 



I , l-(l-xi)X x 



— (V2 



I — X3 



d V2 In f 1 



xi X 



x2(l-x 2 )P 



i = xi x(x2 x 3 ) +x 2 (1 - x 2 ) /3(1 - (1 - x 3 ) y 2l 1 - x{). 
Now we integrate by parts, obtaining: 



(272) 
(273) 

(274) 



1 — Xl 

xiX 

(x + 



I11 1 



xiX 



x 2 (1 - x 2 ) P 



In 1 



,T 3 - X^ 



xi I (1 - a; 3 )x 



In 1 



x 2 (1 - x 2 )/3 2 



x 2 (1 - X2) Pi 



In 1 



six 



x 2 (1 - x 2 ) /3 5 



To compute J s i we perform the change of variable y 2 —* (1 — 2/ 2 )/(l — x 3 ), i.e. 



1 ,1 
^3 / dy 2 

1 J xi 



1 

<k)2 / ^ 3 



1 „l-(l-y 2 )/ Xl 

dy 2 I dx 3 

l-X! JO 



1 - x 3 



(275) 



(276) 
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and obtain: 



Jsi = / dC 3 {{x}) 



1 — (1— x 3 )xi 



dy 2 



For £ s we have the following BST relation: 



£ s =xix + x 2 (l-x 2 )/3 s , f3 s =/3(y 2 ,l 



X 

x *~~~) d x 3 +{yi-{^--xx)X^j d y2 



lnf s >, p = xi+x 2 (l-x 2 ) {l~xi) 2 



Integration by parts gives: 



J si 



2 5, 



dddx} , y 2 ) 



1 — X\ 



1 



1 - x 3 

1 — Xi + Xi 







X 2 / 


Soi 






x 2 J 


£21 J 



-JT X ln|l 
+ (l-a;i) 



1 



In Co 



1 - x 3 
X 



f 2 l ^2 £ 



11 



* ln^ 



■r-2 



ln£ 31 + 2 



*31 



where the new quadratic forms, £ etc, are given by 



£o = ^iX + x 2 (l-x 2 )f3_ , £ 01 =six(i 2 )+i2 (l-x 2 )0 o > 

£i = »i X + ^2 (1 - E2)/?i, Cn = xix{x 2 ) +x 2 (l- x 2 )Pi, 

§ 2 = xix + x 2 (l-x 2 )f3 2 ,_ £ 21 =X!x(x 2 )+ x 2 (l -x 2 ){3 2 , _ 

£3 = xix(0) + x 2 (1 - x 2 ) f3 3 , £, 31 = X!x(x 2 ) + x 2 (l - x 2 ) (3 3 , f3 3 = (1 - xi) 2 x( x 3 ). 

For J21 we perform the same change of variable made for J s \ fEa. (|276|l ). obtaining 



J21 = / dxi dx 2 
Jo 

Secondly, we use Eq. l(5fijl 

x 2 (1 - x 2 ) _ 1 



1 ry2 

dy 2 / dx 3 
Jo 



1 , 1 X 2 (l-^)(2-X!) 

dj/2 / ^3 == • 

l-ai JO J X?s 



x£ s 



5, 



0:2(1-2:2) 1 

+ 3=- I a; 3 



2/3. 



2:2 



X2 (1 - x 2 )(3 6 
xiX 



Similarly for J22 we use: 

2:2 (1-2:2) 1 



x 2 1 



+ i(, 2 - ^) X2 In f 1 + 



2/? s 



3'3 



After integration by parts, we collect the results obtaining: 

1 



J 2 = J 2 \ + J 2 



2B„ 



dCidx}, y 2 ) {2- Xl )I 2 + 



2 5, 



2:1 X 



-/s2, 



/ 2 = -^ln(l 



2:2 (1 - x 2 ) O 



) + 



2:1 X 



/? V ' HX ^ V 2:2 / /?! 

( X3 _ £x ) fi hl f! + x 2 {l-x 2 )-p 2 ,X^l__x 1 ^ / + x 2 (1 - s 2 )%,x 
\ x 2 J (3 2 \ xi \ ) xo B„ V x-\ YlO) / 



2:2 /3 3 
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J. 



s2 



1 x 2 X x + (l-a; 2 )X x +.T 2 (l-2x 2 )(l-X3) 



dC 3 ({a;}) / dy 2 (2 - an) 



'xi xix + x 2 (l - x 2 ) (3 

To compute J s2 we use the same techniques already used for J s \ and obtain 



1 



2 5, 



dC 4 ({x}, y 2 ) (2 -an) 



x 2 X + (l-x 2 )X 



Jsl H Js 



X y ln^ 



I s2 = 2.x 2 (1 - xx) (1 - xz) ln^o + 1 - (1 - sci) x 2 

- x 2 (l-x 1 )X x +x 1 X x ln| 2 -(l-^i)X x ln| 3 + 2x 2 (l-Xi) (1 - x 3 ) 
Summarysing, we write the result for y H as 

(z3=- (vV r ( 1+£ ) (-j-i+w) 



(287) 



(288) 



(289) 



, u2 , dC 4 ({x},y 2 )I° x I%, 



V H = / dC 4 ({x} , y 2 ) 



2B, 



J12 + (2-afi) J a 



4B 2 

x 



9 V a; # 1 x y J 



o (2-x 1 )[x 2 X x + (l-x 2 )X x ]-2 rl l-2x 2 



P 



P 



(290) 



5.8 Collinear limits of V H 

In the previous section we have derived a suitable integral representation for the V H diagram in a generic 
infrared configuration and for arbitrary value of the masses. The resulting representation is stable and can 
be computed numerically even if it contains several polynomial denominators depending on the integration 
variables. This is made possible because each denominator is multiplied by a logarithm or a polylogarithm 
which vanishes exactly at the zeros of the denominator. However it can happen that, for some values of 
the masses, this compensation is delayed and there we encounter numerical instabilities, revealing a region 
where the integrand has strong peaks. This is always the case when one of the two masses in the diagram 
(or both) is vanishing, i.e. in the collinear region. 

For example, if we consider the expression of I® of Eq.J2MJ the second term is proportional to: 



^ x(0) 



— In M 

X3 



-P 2 x 3 + P 2 - M 2 
M 2 



(291) 



In this term the stability around X3 = is at stake when M 2 is small compared to \P 2 \- 

In all these cases one has to compute explicitly the leading behaviour of the diagram, leaving a stable 
remainder. We have considered the following collinear limits of V H : a) m 2 = M 2 <C \P 2 \ and 6) m 2 <C M 2 <C 
|P 2 |. In both cases, we have found that it is better to change Feynman parametrization and proceed in the 
following way (see also section 10.2 ). Starting from Ea. 1)244(1 . we first combine propagators [1] H — [2] H 
with parameter z\, propagators [3]^ — [4]^ with z 2 and [5] H — [6] H with z 3 : 



^ v ° = ^ I dC3i{z}) I ***** [imkmi - 



(292) 



[12] H = [1} H (1 - 21) + [2} H z u [34]* = [4] H (1 - z 2 ) + [3]„ z 2 , [56} H = [5} H (1 - z 3 ) + [6] H z 3 . (293) 
Next we combine [12] H — [34] H with parameter x and integrate in q±, obtaining: 

pi 



7T V 



2e z7r- £ / 2 r 



(2 + j dx j dC^izVixil-x)}- 1 ^/ 2 J d n q 2 



1 



[1234]£ 



2+e/2r 



(294) 
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[1234]„= [q 2 + (l- zi)P2- z 2 P] 



2 . x{zi]pl;m\,ml) x{z2]P 2 ;ml,ml) 



(295) 



x 1 — X 

with x defined in Eq.(|5J|. Finally we introduce the parameter y for the remaining two denominators and 
integrate in q 2 : 



Vo H = ~ 0£) r(2 + e) j\x J\ J dC 3 ({z})[xy(l-x)Y+^(l-y)U- 



2-e 



(296) 



U H = xy (1 - x) (1 - y) [- P 2 (z 2 -z 3 ) (l — zi — zz) + p\ (z 2 ~z 3 ) (1— Z\ — z 3 ) + pj (l-zi-z 3 ) (1- 21—22)] 
+ y(l-x) x(zi\p 2 2 ;mlm 2 2 ) + xy X (z 2 ;P 2 ;ml,ml) + x (1 - x) (1 - y) x(z S ;ph m\, mjj) (297) 

Inserting the values corresponding to the infrared configuration of Fig. lT2l we get (after the change of variable 
Z 3 -> 1 - z 3 ): 

( fx 2 V T(2 + e) 



V H — - 



\P 2 TT 



pi 



dC 5 (x,y,{z}) [xy(l- x)] 1+e/2 (1 - y) (a + a TO e m + a M e M ) 



(298) 



a = xy (1-y) \z\z 2 z 3 + (l-zi)(l-z 2 )(l-z 3 )] + (x + y - xy) z 2 (l-z 2 )| (299) 

a m = x^y (l-x) (l-y) \z x z 2 + z 2 (l-z 3 ) + (l-zi)(l-z 3 )] + (l-a;)(l-y) 2 (l-Z3) 2 + y (x + y - xy) z 2 j 

a M =y\x (l-x) (1-y) [z x z 3 + zi(l-z 2 ) + (l-z 2 )(l-z 3 )] + (l-x)(l-x + xy) z\ + x (x + y - xy) (l-z 2 )\ 

We have also introduced the (small) ratios e m = m 2 / P 2 and e M = M 2 /P 2 . The polynomials <xq, (X\ and a 2 
are positive definite in the integration region, but vanish at the hedge. These zeros are responsible for the 
infrared pole and the collinear divergencies. 

At first we perform a Mcllin-Barnes splitting, once for case a) and twice for case b): 



a) V H 
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b) 



v H — 



P 2 ir 



T(2 + e) / 1 



pi 



1 \ r+ ico 
) / ds dt B (s, 2 + e - s) B (t, s - t) 

2W J -too 



x / dC 5 (x,y, {z}) [xy (1 - x)] 1+t/2 (1 - y) a *(a m e m ) s 2 e (a M e M ) 
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where 0<Res<2 + e and < Ret < Res. For the s and t integrations we close the integration contour 
in the positive real half-plane and compute the residues at the poles of s and t (which are t ~ s ~ n + ae, 
n > 2). Since we are interested in the leading behaviour in e TO — * and e M — > 0, we restrict our attention to 
the poles at n = 2 (higher values of n would lead to terms proportional to e"~ 2 or e^ -2 ). 

A problem connected to case b) is the appearance of poles at t = 2 + e together with poles at t = s 
and s = 2 + e. The residue of these poles in £ generates terms containing all powers of e m /e M ; without an 
assumption on the value of the ratio e m /e M , we would be forced to consider all poles at s = n + ae (n > 2) 
and to resum the series. For this reason we limit our analysis to m 2 <C M 2 and neglect in this way all 
contributions proportional to m 2 /M 2 . 

The poles in s and t come from the integration over Feynman parameters and are, obviously, related 
to the zeros of ao, a,\ and a 2 . These three polynomials can vanish only at the end point of the integration 
region; in order to bring all the zeros to the origin, we first split the integration domain of each variable (for 
x and y this is actually not necessary) according to the following rule: 



dzf(z) = 



1/2 



dz 



dz 



1/2 



/(*) = 



1 



dz 



'(5)+' 



(301) 
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In this way we disentangle z = from z = 1; then we remap each sector into [0, 1], moving all end-point 
singularities in z = 0. Applying this decomposition to the Zi integrals we obtain eight new integrals. Consider 
one of the integrals generated in case a) where the poles are at s = 2 + a e, 



I(s) = / dC 5 (x,y,{z}) (xy) 



l-s+e/2 



(l~x) 1+ ^(l-y) A-° (B al e m y 



(302) 



A al = (1-x) (l-y) [z lZ2 z 3 + (2-2 1 )(2-z 2 )(2-2 3 )] + 2{x + y-xy)z 2 (2-z 2 ) 
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B al = (l-x) \2xy(l-y)(A-2z 3 + z 1 z 3 )+y(l-x + xy)z 2 1 +x(l-y) 2 (2-z 3 ) 2 \ + Axy(x + y-xy) 
To extract the proper behaviour around x = y = 0, we apply a sector decomposition (see section \'6.2\ : 



I( S ) 



dC b (x,y,{z}) (xy) 1 s (l-xy) 



+ x*/ 2 (l-xyy/ 2 (l- y )A-*B al „ 
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A aX ,x = (1-x) (l-xy) [z x z 2 z 3 + {2-z 1 ){2-z 2 ){2-z 3 )\ + 2x(l + y - xy) z 2 (2-z 2 ) 
B alA = 2xy(l-x)(l-xy)(4-2z 3 + z- i z 3 ) + y (l-x)(l-x + x 2 y) z\ 

+ (l-x)(l-xy) 2 (2-z 3 ) 2 + 4x 2 y (1 + y-xy) 
A aX ,2 = (l-xy) (l-y) [z x z 2 z 3 + (2- z{)(2- z 2 ){2- z 3 )] +2y(l + x- xy) z 2 (2-z 2 ) 
B a i,2 = 2xy(l-xy)(l-y)(4-2z 3 + Zl z 3 ) + (l-xy) (l-xy + xy 2 ) z\ 
+x (1- xy)(l-y) 2 (2- z 3 ) 2 + 4xy 2 (1 + x-xy) 
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The poles are at s = 2 + e/2 and s = 2 (both -B i,i and B a \ 2 are not vanishing at x = y = 0). and 
we use Eq.© to extract the poles. All integrals that appear in the evaluation of the diagram have been 
analyzed according to this strategy. The coefficients of the collinear logarithms are computed with a numerical 
integration with a result that agrees with all analytical expansions which have been presented in the literature. 
Our numerical findings can be summarized by means of the following formulas: 



For case a) we obtain 
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5.9 A detailed study of the V H configurations 

We use V" as a prototype for discussing the classification of configurations that, for each family of 
diagrams, are infrared divergent (at least the QED-like). The tools are the corresponding Landau equations 
and a necessary condition for some configuration to be infrared divergent is that the Landau equations are 
satisfied. At first we enumerate internal masses and internal vertices of V H , as done in Fig. Then we 



P2 




Figure 13: Enumeration of internal masses and of internal vertices for V H . 

look for a non-trivial solution of the Landau equations where P 2 or p\{p 2 ) are unconstrained and one or 
two internal masses are zero, according to the rules of the standard model: there is no tri-linear vertex with 
two photon lines. In this operation particular care must be devoted in recognizing topologically equivalent 
diagrams. An example is illustrated in Fig. 1141 where it would be enough to compute our standard definition 
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Figure 14: Equivalence of two configurations of the V H -type containing two internal photonic lines. 

for V H and to perform a permutation of its arguments. The equivalence can be shown by simply enumerating 
the internal vertices and corresponds to (12345) <-> (34521). A second example is given in Fig. 1151 where 
the equivalence is (12345) <-> (14325). For all diagrams but V" we can solve the set of Landau equations, 
looking for the leading or the sub-leading singularities, in the general case and we have an explicit condition 
on internal masses and external momenta which express the singular configurations. For V" there is no 
known general expression and the set of eight equations that we can write must be examined configuration 
by configuration. For instance, given propagators [i] H , i = 1 • • • 6, the standard configuration with p\ = 
—M 2 , P2 = —to 2 , mi = ?Ti3 = to, and TO4 = TO5 — M, m-i = m§ = 0, has a solution with P 2 unconstrained, 
&2,ctQ ^ and ai = = — — 0. Clearly, this has to be put in correspondence with a non-leading 
singularity; we stress that to have a non-trivial solution is the necessary condition for a singularity, but the 
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Figure 15: Equivalence of two configurations of the V^-type containing two internal photonic lines. 

final presence of the singularity follows from infrared power-counting. Indeed, this configuration of V H has 
a single pole at e = 0, as well-known. 

The configuration where mi = — needs not to be considered, given the equivalence of Fig. 1151 This 
exhausts the configurations where we have two internal photons, no tri-linear couplings with two photonic 
legs and both photons attached to a current with one external, on-shcll, line. If the last condition is relaxed 
we can consider a configuration mi = mi = 0, 7713 = m% — m, mi — m§ = M where we only require that 
both lines attached to a photon have the same mass. From the system of Landau equations we observe 
a possible solution with only ai ^ which requires p\ = —M 2 (on-shell condition) and some peculiar 
relation between P 2 and p 2 , namely p\ = —2 (—Mm + 1/2 M 2 — 1/2 P 2 ). Finally, we consider the case 
of only one internal photon. The typical configuration has p\ — —M 2 , p\ — — m 2 , and mi = UI3 = m, 
mi — = Mi, rri2 = 0, rriQ = M2. This is a non-trivial solution but infrared power counting shows that it 
does not correspond to a singularity. One starts from some parametrization of the diagram, performs sector 
decomposition and obtains expressions of the form given in 

Eq.j2Zl with A>0. 

6 Behaviour of two - loop vertices in the collinear limit 

All results derived in this paper have a common property: diagrams are given in terms of integral 
representations where integrands are smooth functions thus allowing for a stable numerical integration. 
These integral representations have been obtained in several ways: we have used BST functional relations 
which force the appearance of B x factors in the denominator; their zeros correspond to anomalous (normal) 
thresholds of the diagram. Far from threshold we register very good stability in the numerical integration, 
but even near threshold we have been able to find quite stable results. 

For those cases where we have not been able to find the proper BST algorithm we, nevertheless, succeeded 
in writing an integral representation of the following form: 

J dC k ({x})j In or J dC k ({x})jU n (^j (312) 

where A, B are multivariate polynomials in the Feynman parameters. These representations generalize the 
class of Nielsen polylogarithms where we only deal with monomials in one variable. In a word, two - loop 
diagrams are always reducible to combinations of integrals of the type given in Ea. (|312|) where the usual 
monomials that appear in the integral representation of Nielsen - Goncharov generalized polylogarithms are 
replaced by multivariate polynomials of arbitrary degree. We have made no attempt towards an analytical 
classification of these new transcendental functions; rather, we compute them numerically, possibly after the 
elimination of apparent singularities by means of sector decomposition techniques. 

Problems with numerical integration arise whenever A and B become small in the same region of Feynman 
parameters and this usually occurs when some mass is much smaller then the largest scale of the diagram, 
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i.e. in the collinear regions. In order to have a better understanding of the collinear limit, we consider the 
case of Eq. l|243[l ; here one of the terms in our representation is 
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B = M 2 y. (313) 



Here numerical instabilities may arise from two specific regions: M 2 is much smaller than A which means 
that B is always small, m 2 is much smaller than |P 2 | so that one zero of A is near y — 0, where also B 
vanishes. Both regions are indeed of collinear nature. To conclude this section we list some useful expansions 
of one-dimensional integrals which often occur in our results. The following formulae are used to extract 
collinear logarithms in order to obtain numerical stability in the collinear regions. 

The expansions have been implemented in the procedure of numerical integration which decides if we 
need to compute the diagram in the collinear region or not and uses the collinear expansion or the initial 
integral representation. 

All collinear divergencies occur at the border of the integration domain which, for one-dimensional inte- 
grals means x = or x = 1. The latter case can always be eliminated with a change of variable, x — > 1 — x. 
Therefore, all integrals in collinear regions which occur in our results can be cast in one of the following 
forms, where Q(x) — ax 2 + bx + c — iS, with £ ^ 0+ and R(x) = Q(x)/(dx): 
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dx 
o Q{x) 



\nR{x) \n{l-R(x)). 



In our example we find a /| and the two unstable cases correspond respectively to d — > and c — * 0. 

When d — ► 0, the behaviour of 7q is trivial and I\ is regular. Furthermore, to extract the divergency for 
l\ and l\ it is enough to use well-known relations among (poly)logarithms of argument z and 1/z and then 
apply the BST method fEa.l|S7|l and Eq.lJSJl). After integration by parts, the logarithmic behaviour in d 
arises naturally without any expansion. The results, collected in appendix^] are therefore valid everywhere. 
In order to extract the divergent behaviour for c — * 0, we proceed in the following way: 



dxf(x,Q(x))= / dx[f(x,Q(x)) — f(x,bx)]+ / dx f(x, bx + c) + / dx [f(x, bx) — f(x, bx + c)], (315) 



where / is one of the (poly)logarithms of Eq. l|315|l . Only the second integral in Eq. l|315|l is divergent when 
c — * 0, we can set c = in the first (subtracted) term and the third one is O (c). The results for the expansion 
are again listed in appendix lEl 



7 Numerical Results 

In this section we present numerical results for infrared configurations corresponding to two-loop three- 
point scalar functions. All results are computed numerically using the analytical expressions given in this 
paper. A particular attention has been devoted to collinear regions where the mass of the outgoing particles 
is small compared to the incoming momentum. 

In our general program, aimed to a numerical evaluation of multi-loop, multi-leg Feynman diagrams we 
have developed a set of FORTRAN/77 codes which go from standard Aq, ... ,Do functions to diagrams 
presented in this paper. This new ensemble of programs which includes the treatment of complex poles 32 
will succeed to the corresponding Library of T0PAZ0 [231 ■ 

The whole collection of codes is heavily based on the NAG- library [5D|: while completing the analysis it 
became natural to consider a migration of the whole set of programs to a stand-alone, FORTRAN/95 version. 
The outcome of our decision is a brand-new version of our numerical code, CoopBack |34| . which is fully 
based on quasi - Montecarlo methods (with a possible extension to a parallelized version) and is presently 
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under construction. Our FORTRAN/95 version has been adapted from the automatic multi-dimensional 
integration subroutine DKBVRC, written by Alan Genz; DKBVRC uses randomized Korobov rules |2E] 

In the migration to FORTRAN/95 we have experienced huge gains in CPU-time, a better numerical 
precision and we foresee a future release of the code, although the numbers produced for this paper still rely 
on the (by now) old version. For present numerical results all the vertices are evaluated using the routine 
D01GDF or D01EAF j^O]- The first routine is based on the Korobov-Conroy [2H| number theoretic method 
with a Monte-Carlo error estimate, while the second one uses an adaptive subdivision strategy. 

For those configurations where we could compare with analytical results we have adapted our setup to 
match the known examples. Alternatively we have considered few physically relevant cases, selected among 
those not presented in the literature, and extracted from processes like Z* — > //, H* — > W + W~ , t — » W + b, 
etc. In this paper we use the following input parameter set: 

M w = 80.380 GeV, M z = 91.1875 GeV, M H = 150 GeV 
m t = 174.3 GeV, m 6 = 4GeV, m e = 0.510999 MeV (316) 

Comparisons have been performed with (T^j and the results are shown in Tab. Tab. EH Tab. 03 Tab. |H1 
Tab. 1101 Tab. 1121 and Tab. 1131 with ^2] an d the results are shown in Tab. 1141 with ^3] an d the results are 
shown in Tab. 1151 In the comparison with the work of Bonciani - Mastrolia - Remiddi ^5], the evaluated 
finite part is defined through the following e expansion 2 : 

Vav = V-f \ + V~j - + V° T , etc. (317) 

For the comparison with the work of Davydychev - Kalmykov |12j we use V* of Fig. 1101 with mi = 7772 = 
7773 (= 180 GeV), rrii — 0, i > 4 and pf = (which is related to H* — > gg). We compute 

V a = (^) T 2 (l + I) (V-2;k i> + V-i )X \ + V . K ). (318) 

Note that the presence of a double pole in Eq. l|318|l reflects a collinear singularity besides the infrared one. 
For the comparison with the work of Davydychev-Smirnov ^2] we consider the diagram V* of Fig. 1101 with 
m\ = 777,4 = m an d m 2 = m 6 = M. The analytical calculation is valid in the region m 2 <C M 2 , \ P\ 2 , and we 
define 

V C K = (—) r (2 + e) (V-2-.k ^ + V-X;k ~ + v o-,k)- (319) 

For all comparisons performed in this paper we have found an excellent agreement with analytical calcula- 
tions, 3 therefore signaling a satisfactory status of the overall goodness of our numerical algorithms. Needless 
to say, we have been able to produce results having no counter-examples in the literature. 

In Tab. 1161 we give a sample of results for the numerical integration of tensor integrals. The relevant 
message is that all analytical expressions which have been implemented in our code for tensor integrals has 
been derived using the same techniques already used for scalar configurations. 

The topology chosen for tensor integrals is V c - K , which enters in the computation of the fermionic cor- 
rections to sin 2 o ff and these results have been already used in ; this emphasizes the overall relevance of 
our results. 

The tensor structure considered in this brief example is: 

_ M 2e f d"gid"g 2 y, _M 2£ [ dn , gf _y y , u 

// 2e r n v 2 



d^ = ^4 / dnqidnq2 W 2 W= ^I^IW.J^/i+W*' 8 '. ( 320 ) 

11/— l i \ K ii,j2=i 

where [l] K has been defined in Ea. (|148l) . The e expansion for the coefficients is V c - K — /e + V° K . 

2 Note that these authors define n = 4 — 2 e. Therefore the e of the following formula corresponds to 2 e of 1151 

3 In several cases we had to code the analytical results since the authors did not present explicit numerical results 



48 



8 Conclusions 



In this paper we have analyzed a special component of our project aimed to a numerical evaluation 
of physical observables at the two-loop level: since QED/QCD are an integral part of any realistic theory 
and they are characterized by the exchange/emission of massless gauge-bosons leading to infrared divergent 
parts we had to prove that infrared and also collinear configurations can be treated within the same class of 
algorithms which we have used for massive configurations or, at least, within some simple extension of them. 

The main result, therefore, has been to assemble relatively simple expressions for scalar two- loop vertices 
in a systematic and coherent manner so that they can be used for practical calculations. Our results introduce 
integral representations which are well suited for numerical integration and represent a generalization of the 
familiar Nielsen - Goncharov multivariate poly logarithms. 

Confining most of the paper to scalar configurations should not be confused with a limitation of the 
method. Tensor integrals that arise in any non-trivial theory, due to the spin structure, simply add extra 
polynomials of Feynman parameters in the integral representations of the diagrams. One can easily see, for 
instance from Ea. (|73|) . that these additional polynomials can only change the numerator structure inside 
Ea. i|7tjf) which is at the basis of our results: we simply get an hypergeometric function with a different list of 
arguments and the whole derivation can be carried through along the same lines used for a scalar integrand. 

Another important issue that has been addressed in this paper concerns the systematization of any 
procedure for implementing infrared divergent graphs in a realistic calculation, at least when using a modern 
language and when considering QED (and also QCD) as embedded into a larger theory. Indeed, QED alone 
with a massive regulator has been treated long ago in a seminal paper by Cvitanovic and Kinoshita |24|. 

In the final stage of our project we will generate diagrams with the help of QraphShot each dia- 
gram where an internal photon (or gluon) line appears will be subjected to a special investigation, namely 
the corresponding Landau equations will be examined. As soon as they are fulfilled by the kinematical 
configuration that we are scanning and as soon as the filter of infrared power-counting is passed we know 
that the configuration is infrared singular and the appropriate subroutine will be initialized returning nu- 
merical answers for the residues and the finite part. Understanding these motivations will hopefully explain 
our preference for extending the numerical treatment to infrared divergent configuration despite the recent, 
spectacular, progress in analytical evaluation. 

For instance, one of our configurations, V K , has been already computed by Davydychev and Smirnov 
in the limit m 2 <C M 2 , \P |; we have found excellent agreement and have able to extend the numerical results 
to all values of m, therefore allowing for QCD corrections to the top quark decay without approximations for 
internal masses. For another setup of the same V K configuration we have another result by Davydychev and 
Kalmykov |12j which is relevant for Higgs decay into two photons or two gluons; again we found excellent 
agreement and are able to produce numerical results for the same diagram embedded into the standard model, 
e.g. also with VK-lines and not only gluons outside the inner massive quark loop. Furthermore, we have 
been able to perform a numerical test of several analytical results by Bonciani, Mastrolia and Remiddi 15 
corresponding to the whole set of two-loop topologies. 

To summarize, we have been able to present all formulas that form the basis for numerical evaluation 
of infrared residues and infrared finite parts of arbitrary infrared configurations of two-loop vertices. The 
language may sound unfamiliar but our results have far reaching consequences; for instance, some of the 
results presented here have already been used for computing two - loop electroweak pseudo - observables 7 . 
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A Taylor and Laurent expansion of Euler's functions 



Here we collect results needed in expanding Euler's functions; 7 is the Euler constant and £(x) is the 
Riemann zeta function. 



r(l + z) = l + jr G n z n , r(z) = T{1 + z) =- + J2G n 



n=l 



n=0 



Gi=-% G 2 = ^ 7 2 +C(2) 



G3 = - 6 



7 3 + 3C(2) 7 + 2C(3) 



G 4 = ^ 7 4 + ~ C(2) 7 2 + ~ C 2 (2) + I C(3) l + \ C(4). 
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B Nielsen polylogarithms 

Throughout the paper we have used (n , p are positive integers) 



(n - 1) 



>-l r 1 

— - / dx — In™ - x ln p (l — zx), S n -i \{z) = Li„ 
P' Jo x 
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C Properties of the hypergeometric function 

The Gauss hypergeometric function |27] is defined by: 



iFi(a, b; c; x) 



r( c ) 



dzz 6 - 1 (l-z) c -"- i (l-a;z)- a 



Re & > Re c> 0. 



(323) 



r(6) r( c -6) y 

The special case c = b + 1 frequently occurs in this paper. Sometimes the Gauss hypergeometric series is 
used (circle of convergence |x| = 1): 

T(a + n) b 



2*1(0, b; b+l; x) =Y^ 



n=0 



r (a) & + n T (n + 1) ' 



(324) 



An important property of 2-F1 used in our paper is (1 — c, b — a, c — b — a are not integers): 

3 F 1 (a , b; b + 1 ; x) = -A- (-x)" 2 F 1 (a, a - 6; a - b + 1 ; -) + r ( b + 1 )^( a ^ b ) ( _ a;) -6 ) 
— a x 1 (a) 

where |arg(— x)| < tt. In other cases we need an expansion around a vanishing e; few examples are listed 
below (see also [35]): 
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In several cases the hypergeometric function (with integer arguments) leads to elementary transcendental 
functions: 



2*1(1, 2; 3; z) 



- ln(l - z) + 1 



2 Fi(2, 2; 3; z) = - { - ln(l - z) + 1 



1 



1 



2^(1, 1; 2; z) = -- ln(l-z), 



2 Fi(2, 1; 2; z) = 



1 - z 



(327) 
(328) 



D Computation of the integral J^(i) 
The function J^(i) is defined by: 

In a, (a;) In 71 aj(t/) 
ai(ar) a^y) 



C(0= / / dyln n xln m y 



ai(z)=x(z), a. 2 {z)^x( z ), (329) 



where \ and % arc defined in Eq. and Eq. Hll|) , respectively. The integral in Eq. (|329() can be computed by 
using BST functional relations. We first use the following relation: 



aT\y) \n k a t (y) = ±- [ln fe a t (y) - d y (ln fc+1 a,(y) - ^+^(0)) 



(330) 



where Zi = X and Z 2 = X (see, once again, Eq.©). With respect to the usual BST relation (see Eq. 1151(1 ) 



with have added a term, d y In cti(0), which is actually zero. After the integration by parts, we treat the 
x-integral through another BST relation: 

a-\x) hx k ai (x) = -j- \ln k ai (x) - *~f* d x \n k+1 a^x) 

y Z(K + 1) 



(x - Zi)a i 1 (x) ln fc oti{x) 



The second integration by parts gives: 



1 a 



x i fc+1 



2P 2 fc+1 



In ai(x). 



J, 



(«) = 



45^ 



fc-1 , 



R 

nm 



(0 , ^ 



(0 



P 2 (fc + ft) 



(331) 



(332) 



where, new, auxiliary quantities have been introduced: 



= jf <ic J dyL^iy) In"" 1 x | In a; /^(x) [Z,(x) + 2fc] + n (x - Z ( ) 

- X! ! r(- 1 ) m " PZ ' i (°) / ^ ^ln" +p x^- 1 (x)[/ l (a:) + 2fc] 
p=0 p ' ^° ^ 

+ (In* + n + p) (x-Z t ) ln'^- 1 xZ*(x)j - JdyL^iy) 

6 n , (1 - Zi) 11(1) - (y - Zi) In" yl\ (y)] + S n<0 (-1)" 1 ml (1 - Z-) l*(l) l*(0) 



fc-l , h-l 



(i) = (n + to) / In 



n+m — 1 



^+' l (x) 



Furthermore, we have defined: 



i* +h (i)-(i + l)ii(i)iHo) + li? +h (o) 



l i (x)=hxa i (x), L k -\y) = \n n yl h i - l {y) [k(y)+2k] +n(y-Z i ) In™" 1 y 



ijKv) 



(333) 



(334) 
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E Useful expansions 



Here we discuss various expansions of I ' n (denned in Section [fJJ) for c — ► in terms of polylogarithms 
(Eq.J222J0- Let Q(x) = ax 2 + bx + c and e = a + b + c; we have 

i — h 



Io' n (a,b, c, d) 



In" 



h=l 



1=0 



l\ (n-h-l)\ 



kl 



(n - h)\ 



In" 



- ] Sk+:i /, 



ln n-h-l 



0(c). 



Then we consider a class of integrals which generalizes Eq. (|130|) : 



Jn(a,b, c,d) 



dx 



In" x k Q(x) 



Q(x) d ' 

They can be computed by using the BST relation Ea. (|51|l which in the present case reads as follows: 



Q{x) 



In 



c Q(x) 



1 

B 



In' 



Q(x) 



x-X 
2{k + 1) 



cL In 



fe+ 



i Q(x) 



X ~~2^' 



fe 2 -4, 
4T" 



(335) 



(336) 



(337) 



After integration by parts, we get: 



J2 = I f dx (n In" 1 x + In™ x) In 

2B {J 



Q(x) 



nXI^ 1 ' 1 + 2(-l) n r(n+l) 



jk-l _ 



1 



rfc-1 



2B k 
1 

2Bk 



dx 



2 k ln n x In* -1 ^ + (n In™" 1 * + ln n x) ln fe ^ 



dx In 



in 



d 

Q(x) 



-nXL 



n — l.k 



2k 



(1 - X) ln fc - - X ln fc 



where for the first two results n > 1. Now we give the expansion of 7™ for c — » 0: 



'^Ind 

x 



1 dx 



LL, 



Q(x 
dx 



= -Lia - 



6 + d 



Li 2 -- 



In - In 1 + - 



(338) 



•0(c), 



1 dx 



■ Li. 



ax + b 



In 



LL 



S n -U 



0(c), (339) 



with n> 2. Next we give the results for an d f° r d ^ 0. The final expression has been derived following 
the strategy sketched in Section El and is valid also for d far from zero. 



4 



T 2 



dx 

o Q( x ) 



In 1 



Q{x) 
dx 



(l-X)Li 2 ( -- ) +XI\ 



-2B 1 ldX 



T 1 -4- T° 



In 1 



Q(x) 



Li 2 - 



dx 
Q(x) 



dx 
o Q( x ) 



■ Lio 



Q(x) 



l 

2B 



dx 



Lio 



dx 
Q{x) 



Li 3 



dx 
Q(x) 



{ i-x)u 3 -^)+xi!\-a2)j° -^Ji- 1 -Jl + Jl 



(340) 



Finally, expanding d^ 1 for c ^ we obtain: 



lo^r -Li 2 - 



6 + d 



1 dx 



o 

In — Li,, 

c 



Li r , 
& 

"d 



-Li a I - 

ax + 6 V 



d 



d 

Li n +i 



d 
c 

a + 6 



In - In ( 1 + - 



+ (n + 1) Li n+ i 



iln^l 



0(c) 



- -2Li 2 -- 



O(c) , 



n > 2 



(341) 
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F Tables of numerical results 

In this section we collect a sample of our numerical results, based on the setup of Eg. l|316[) . 





s 


Re V£. t 


Ini V 

AAA1 v a;i 


(Jur 


[500 GcVJ 


— 2.4142(1) X 10 


4.12(1) x 10 


BMR 




-2.414166 x 10" 2 


4.121472 x 10~ 3 


Our 


Ml 


-0.61103(4) 


9.447(6) x to" 2 


BMR 




-0.6110145 


9.446769 x 10~ 2 


Our 


100 m 2 


-6.04493(6) x 10 7 


3.0986(2) x 10 7 


BMR 




-6.044903 x 10 7 


3.098533 x 10 7 


Our 


4.0001 m 2 e 


-6.878781 x 10 10 


3.906172 x 10 11 


BMR 




-6.878781 x 10 10 


3.906172 x 10 11 


Our 


ml 


8.533838 x 10 8 





BMR 




not available 





Our 


-TOO to 2 


5.463023 x 10 7 





BMR 




5.463023 x 10 7 






Table 1: Diagram V a; i. Comparison with the results of I15| (BMR). The setup for Fig.|^|is: mi = mi — m = M = 
m e . The unit of mass fi is 1 GeV. Only the infrared finite part is shown. The results are in GeV -2 . Unless indicated 
our relative error is below 10 -7 . 



Process 


s 


M, m 


mi, m2 




v a;l 


Z* -> e+e" 


[500 GeV] 2 


m e , m e 


mb, mb 


-4.7574(4) x 10- 2 


2.55694(6) x lO" 2 




[500 GeV] 2 


m e , m e 


mt, m t 


-1.147833(1) x 10" 2 


4.303633(8) x lO" 3 


Z -> e+e- 


Ml 


m e , m e 


m b , m b 


-7.91577(3) x 10- 2 


4.74278(5) x 10~ 2 




M 2 

Z 


m e , m e 


m t , m t 


-0.2243587 


9.752214(1) x 10~ 2 


Z* -> e+e- 


100 m 2 


m e , m e 


m b , m b 


-5.665645 x 10 6 


9.745030 x 10 5 




100 m 2 


m e , m e 


mt, m t 


1.987627 x 10 5 


4.155691 x 10 5 


W e+v e 


K 


m e , m v 


mt, m t 


-0.7168993 


0.1227244 


W* -> tb 


[200 GeV] 2 


m t , m b 


m b , m b 


-7.158(2) x 10- 3 


3.3638(2) x 10~ 2 


t -> W+b 


2 


m b , M w 


M w , M w 


-2.709610(3) x 10- 2 


3.182718(2) x 10" 2 



Table 2: Diagram V a; i The setup refers to the diagram of Fig. |S|The neutrino mass has been arbitrarily set to 
m v = 0.1 eV The unit of mass fi is 1 GeV Only the infrared finite part is shown. The results are in GeV -2 Unless 
indicated our relative error is below 10 -7 . 
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5 




im V b-I 


111 


fcnn PnVl 2 
[OUU KjC V J 


7 1 (197R9 v 1 n — 2 
— /.lUZfOZ X 1U 


1 Q1 AOfoO w 1 n — 2 


XJlVLrl 




— / .lU/ooO X 1U 


1 01 /iog^ v 1 n — 2 
i.yi4zyo x iu 


Our 


M 2 






BMR 




-1.590738 


0.4462335 


Our 


100 m 2 


-7.073558 x 10 7 


4.1480(2) x 10 7 


BMR 




-7.073558 x 10 7 


4.148499 x 10 7 


Our 


4.0001 m 2 e 


-7.264768(2) x 10 10 


4.364353 x 10 11 


BMR 




-7.264767 x 10 10 


4.364353 x 10 11 


Our 


mi 


1.000079 x 10 9 





BMR 




not available 





Our 


-100 ml 


6.714222 x 10 7 





BMR 




6.714222 x 10 7 






Table 3: Diagram Vb i. Comparison with the results of |15| (BMR). The setup for the diagram of Fig. is: 
mi = m — M — m e , 7712 = 0. The unit of mass /j, is 1 GeV. Only the infrared finite part is shown. The results are 
in GeV -2 . Unless indicated our relative error is below 10 -7 . 



Process 


s 


M, m 


mi, m 2 




Im V b ° :I 


Z* -> e+e" 


[500 GeV] 2 


m e , m e 


m e , M z 


-8.675940 x 10~ 3 


4.466520 x 10~ 3 




[500 GeV] 2 


m e , m e 


m u , M w 


-8.462338 x 10~ 3 


4.491430 x lO -3 


Z -> e+e" 


Ml 


m e , m e 


m e , M z 


-0.1618863 


7.892145 x 10~ 2 




Ml 


m e , m e 


m v , M w 


-0.1551424 


7.718541 x 10~ 2 


Z* -> e+e" 


100 m 2 


m e , m e 


m e , M z 


-5.936619 x 10 5 


-9.349849 x 10 4 




100 m 2 


m e , m e 


m v , M w 


-7.285072 x 10 5 


-1.167232 x 10 5 


W* -> tb 


[200 GeV] 2 


mt, mb 


m bl M w 


-2.383592 x 10~ 2 


5.899083 x 10~ 2 




[200 GeV] 2 


m t , m b 


m t , M H 


-3.364485 x 10~ 2 


7.477385 x 10~ 2 


H* -> W+W- 


[200 GeV] 2 


M w , M w 


M w , M H 


-7.102816 x 10^3 


3.303713 x 10~ 2 



Table 4: Diagram Vb ; i. The setup refers to the diagram of Fig. The neutrino mass has been arbitrarily set to 
m v = 0.1 eV. The unit of mass /i is 1 GeV. Only the infrared finite part is shown. The results are in GeV" 2 . Unless 
indicated our relative error is below 10 -7 . 
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s 


ReV° 


' a;M 


Our 


[500GeV] 2 


2.0928(5) x 10 3 


-1.95(3) x 10 2 


BMR 




2.092863 x 10 3 


-1.941897 x 10 2 


Our 


Mf 


5.659(1) x 10 4 


-5.84(3) x 10 3 


BMR 




5.658944 x 10 4 


-5.837323 x 10 3 


Our 


100 m 2 


6.28283(2) x 10 12 


-1.70618(2) x 10 12 


BMR 




6.282796 x 10 12 


-1.706180 x 10 12 


Our 


4.0001 m 2 e 


-1.447472 x 10 20 


1.6161(6) x 10 21 


BMR 




-1.447472 x 10 20 


1.616432 x 10 21 


Our 


m 2 e 


-3.11237240 x 10 14 





BMR 




not available 





Our 


-100ml 


-6.544718 x 10 12 





BMR 




-6.544718 x 10 12 






Table 5: Diagram V a;M . Comparison with the results of |15| (BMR). The setup refers to the diagram ofFig. [7| 
mi = m.2 = m — M = m e . The unit of mass fi is 1 GeV. Only the infrared finite part is shown. The results are in 
GeV -4 . Unless indicated our relative error is below 10 -7 . 



Process 


s 


M, m 


mi, m2 




lmV° 

y a;M 


Z* -> e+e- 


[500 GeV] 2 


m e , m e 


m b , m b 


91.02032 


1.5233(1) x 10~ 6 




[500 GeV] 2 


m e , m e 


mt, mt 


-24.61950 


-9.6641(1) x 10" 10 


Z -> e+e" 


Ml 


m e , m e 


m b , m b 


2.736577 x 10 3 


4.4153(1) x 10~ 5 




Ml 


m e , m e 


m t , m t 


-7.401994 x 10 2 


-2.114765 x 10" 6 


Z* -> e+e- 


100 m 2 


m e , m e 


m b , m b 


8.061021 x 10 11 


3.308893 x 10 10 




100 m 2 


m e , m e 


mt, m t 


-2.392478 x 10 11 


-4.086825 x 10 10 


W* -> tb 


[200 GeV] 2 


mt, m b 


m b , m b 


-6.5246(4) x 10" 5 


3.5182(7) x 10- 5 




[200 GeV] 2 


mt, m b 


mt, mt 


-6.132929 x 10~ 5 


-6.549233 x 10~ 6 


H* -> W+W- 


[200 GeV] 2 


M w , M w 


M w , M w 


-1.493808(2) x 10~ 7 


-3.513014(2) x 10" 7 


t -> W+b 


ml 


M w , m b 


m b , m b 


-3.5328(2) x 10- 5 


1.5688(4) x 10~ 5 




2 

m t 


M w , m b 


M w , M w 


-2.341537 x 10~ 5 


-2.486382(2) x 10" 7 



Table 6: Diagram V a; M- The setup refers to the diagram of Fig. □ The unit of mass /i is 1 GeV. Only the infrared 
finite part is shown. The results are in GeV -4 . Unless indicated our relative error is below 10 -7 . 
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Process 


s 


M, m 


mi, m 2 






Z* -> e+e" 


[500 GeV] 2 


m e , m e 


m e , M z 


-25.86641 


2.008353 x 10~ 8 




[500 GeV] 2 


m e , m e 


m v , M w 


-29.73137 


2.488703 x 10~ 8 


Z -> e+e" 


M| 


m e , m e 


m e , M z 


-7.776883 x 10 2 


2.316778 x 10" 6 




Ml 


m e , m e 


m v , M w 


-8.938906 x 10 2 


2.416842 x 10~ 6 


Z* -> e+e- 


100 m 2 


m e , m e 


m e , M z 


-2.171299 x 10 11 


8.580100 x 10 9 




100 m 2 


m e , m e 


m v , M w 


-2.520681 x 10 11 


6.108271 x 10 9 


W* -> 46 


[200 GeV] 2 


m t , m b 


m b , M w 


-2.248708 x 10~ 7 


1.048569 x 10~ 6 




[200 GeV] 2 


m t , m b 


m t , M H 


-4.763299 x 10~ 7 


1.548612 x 10~ 6 


ff* -> Vt^+T^- 


[200 GeV] 2 


M w , M w 


M w , M H 


1.109640 x 10~ 7 


4.589289 x 10~ 7 


i -> VF+6 


mf 


M w , m b 


m e , m v 


3.567051 x 10~ 7 


-2.586087 x 10~ 7 




m\ 


M w , m b 


M w , M z 


-1.412470 x 10~ 7 


3.224342 x 10~ 7 



Table 7: Diagram Vb :M . The setup refers to the diagram of Fig. The neutrino mass has been arbitrarily set to 
m v = 0.1 eV. The unit of mass fj, is 1 GeV. Only the infrared finite part is shown. The results are in GeV -4 . Unless 
indicated our relative error is below 10 -7 . 





s 


Re V° 


ImV° 

iul y c;M 


Our 


[500 GeV] 2 


-1.819980 x 10 5 


4.273610 x 10 4 


BMR 




-1.820011 x 10 5 


4.273685 x 10 4 


Our 


M| 


-4.194142 x 10 6 


1.086325 x 10 6 


BMR 




-4.194142 x 10 6 


1.086325 x 10 6 


Our 


100 ml 


-9.234450 x 10 13 


8.726096 x 10 13 


BMR 




-9.234450 x 10 13 


8.726096 x 10 13 


Our 


4.0001m 2 


1.449091 x 10 20 


-1.754699 x 10 21 


BMR 




1.449091 x 10 20 


-1.754699 x 10 21 


Our 


ml 


1.377881 x 10 15 





BMR 




not available 





Our 


-100 m 2 


1.038162 x 10 14 





BMR 




1.038162 x 10 14 






Table 8: Diagram V C :m. Comparison with the results of [15) (BMR). The setup referring to the diagram of Fig. 
is: m = M = m e . The unit of mass \i is 1 GeV. Only the infrared finite part is shown. The results are in GeV -4 . 
Our relative error is everywhere below 10~ 7 . 
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Process 


s 


M, m 


mi, m 2 , m 3 




imy» K 


Z -> e+e~ 


[500 GeV] 2 


m e , m e 


m e , m e , Af z 


-8.6239(5) x 10~ 9 


4.86523(5) x 10~ 8 




[500 GcV] 2 


m e , m e 


M w , M w m f 


-1.706(4) x 10" 8 


7.109(3) x 10~ 8 


Z -> e+e" 




m e , m e 


m e , m e , M z 


5.11313(2) x 10~ 6 


5.05866(2) x 10~ 6 




M 2 

Z 


m e , m e 


M w , M W! m v 


4.552503(7) x 10~ 6 


-1.269944(6) x 10" 6 


Z* -► e+e" 


100 m 2 


m e , m e 


m e , m e , M z 


-6.700055 x 10 3 


1.226283 x 10 3 




100 m 2 


m e , m e 




-86.21982 


-11.84965 


T^* -> i& 


[200 GeV] 2 


m t , mfc 


mb, m t , M z 


1.0060(2) x 10~ 6 


-6.2131(5) x 10~ 7 




[200 GcV] 2 


m t , m 6 


M w , M H , mt 


1.894547(1) x 10~ 7 


-3.017159 x 10~ 7 


if* _> W+W~ 


[200 GeV] 2 


M w , M w 


m t , m t , m b 


4.216966(3) x 10~ 8 


-1.225831 x 10~ 7 


t -> 


mf 


M w , m b 


M w , m b , m t 


3.4389(8) x 10~ 7 


6.082(2) x 10~ 8 




m| 


M w , m b 


m b , M w , M z 


6.3629(3) x 10~ 7 


3.0576(8) x 10- 7 



Table 9: Diagram V a;K - The setup refers to the diagram of Fig. I1UI The neutrino mass has been arbitrarily set to 
ni v = 0.1 eV. The unit of mass /i is 1 GeV. Only the infrared finite part is shown. The results are in GeV" 4 . Unless 
indicated our relative error is below 1CV 7 . 





s 




ImV^ 


Our 


[500 GeV] 2 


-2.6228(2) x 10" 4 


8.7973(7) x 10~ 5 


BMR 




-2.622920 x 10" 4 


8.797302 x 10~ 5 


Our 


Mf 


-1.64386(8) x 10~ 3 


5.8676(7) x 10~ 4 


BMR 




-1.6438612 x 10~ 3 


5.867706 x 10~ 4 


Our 


100 m 2 


-3.8356(4) x 10 12 


7.7519(4) x 10 12 


BMR 




-3.835560 x 10 12 


7.751602 x 10 12 


Our 


4.0001m 2 


1.844(3) x 10 20 


unstable 


BMR 




1.842500 x 10 20 


5.160917 x 10 19 


Our 


m 2 e 


-1.177792 x 10 15 





BMR 




not available 





Our 


-100 ml 


-5.682692 x 10 12 





BMR 




-5.682692 x 10 12 






Table 10: Diagram Vb ;K . Comparison with the results of [K] (BMR). The setup referring to the diagram of Fig. 1101 
is: m — M = m e . The unit of mass fi is 1 GeV. Only the infrared finite part is shown. The results are in GeV -4 . 
Unless indicated our relative error is below 10 -7 . 
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rlUCcDo 


5 


IK/ , A/i 


mi, TO2 , W3 


KC V C . K 


im V C . K 


ry -f- ^~ 

Zj — > e e 


[oUU CjC VJ 


m e , m e 


m e , M z , m e 


— z.looo^e J X 1U 


0.8 1 (zj X 1U 


r-r _|_ 

Z —> e+e 


A/2 


m e , m e 


m e , M z , m e 


-5.0064(3) x 10- 5 


1.862(2) x 10~ 6 




9 

«e 


M z , m e 


m t , m t , m t 


-2.656960 x 10" 7 


2.424(5) x 10 -20 




9 


M z , m e 


mb, m h , m 6 


6.8986(3) x 10- 6 


-6.3053(4) x 10 -6 


Z* -> e+e" 


100 777 2 


m e , m e 


m e , M z , m e 


-7.19834(2) x 10 3 


3.73760(1) x 10 3 


1U* -> tb 


[200GeV] 2 


m t , m b 


m t , M z , m t 


1.43046(3) x 10- 7 


-2.256417(8) x 10~ 7 




[200GcV] 2 


m t , nib 


M w , m 6 , M w 


5.962(5) x 10~ 7 


6.528(3) x 10- 7 


ff* -> T^+M^- 


[200GeV] 2 


M w , M w 


m t , m b , m t 


3.62055(7) x 10~ 8 


-1.075947(8) x 10~ 7 


t -> TU+fe 


m? 


M w , m b 


m t , 777 b , m 4 


1.22769(4) x 10~ 7 


-1.084614(6) x 10~ 7 




mf 


M W1 m b 


M w , M H , M w 


1.6920(4) x 10~ 7 


-1.49596(9) x 10 -7 



Table 11: Diagram V c;K . The setup refers to the diagram of Fig. 1101 The unit of mass /1 is 1 GeV. Only the 
infrared finite part is shown. The results are in GeV -4 . Unless indicated our relative error is below 10~ 7 . 





s 






Our 


[500 GeV] 2 


1.58634(5) x 10 5 


-3.868(8) x 10 4 


BMR 




1.586339 x 10 5 


-3.878825 x 10 4 


Our 


Ml 


3.61488(4) x TO 6 


-9.75(2) x 10 5 


BMR 




3.614876 x 10 6 


-9.774586 x 10 5 


Our 


100 m 2 


6.7595(2) x 10 13 


-7.0650(4) x 10 13 


BMR 




6.759415 x 10 13 


-7.065317 x 10 13 


Our 


4.0001m 2 


-1.680564(3) x 10 17 


1.1570(6) x 10 17 


BMR 




-1.680565 x 10 17 


1.156966 x 10 17 


Our 


m z e 


-7.970954 x 10 14 





BMR 




not available 





Our 


-100 m 2 


-7.782868 x 10 13 





BMR 




-7.782868 x 10 13 






Table 12: Diagram Vd ;K . Comparison with the results of jT5] (BMR). The setup referring to the diagram of Fig. 1101 
is: m — M = m e . The unit of mass fi is 1 GeV. Only the infrared finite part is shown. The results are in GeV -4 . 
Unless indicated our relative error is below 10~ 7 . 
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s 


Re V£ 


T T rfl 


Our 


500 2 GeV 2 


-1.342815 x 10~ 6 


1.052948 x 10~ 6 


BMR 




-1.342815 x 10" 6 


1.052948 x 10~ 6 


Our 


Mf 


-4.562983 x 10~ 4 


5.478876 x 10~ 4 


BMR 




-4.552983 x 10~ 4 


5.478876 x 10~ 4 


Our 


100 m 2 


2.801(1) x 10 11 


-2.0846(8) x 10 12 


BMR 




2.801721 x 10 11 


-2.084294 x 10 12 


Our 


4.01m 2 e 


-1.87(1) x 10 16 


2.185(32) x 10 16 


BMR 




-1.866708 x 10 16 


2.152803 x 10 16 


Our 


9 

mi 


1.424912 x 10 14 





BMR 




not available 





Our 


-100 ml 


1.21504(6) x 10 12 





BMR 




1.214934 x 10 12 






Table 13: Diagram V H . Comparison with the results of II5| (BMR). The setup referring to the diagram of Fig. 1121 
is: m — M = m e . The unit of mass \i is 1 GeV. Only the infrared finite part is shown. The results are in GeV -4 . 
We have found numerical instabilities for values of s too close to the normal threshold, Unless indicated our relative 
error is below 10 -7 . 





[GeV] 


ReV 0;iC [GeV" 4 ] 


lmV . K [GeV" 4 ] 


Our 
DK 


400 


5.1343(1) x 10~ 8 
5.13445 x 10~ 8 


1.94009(8) x 10" 8 
1.94008 x 10~ 8 


Our 
DK 


300 


5.68801 x 10~ 8 
5.68801 x 10~ 8 


-1.61218 x 10" 8 
-1.61218 x 10~ 8 


Our 
DK 


200 


9.36340 x 10~ 8 
9.36340 x 10~ 8 


-2.84232 x 10~ 8 
-2.84232 x 10~ 8 


Our 
DK 


100 


2.94726 x 10~ 7 
2.94726 x 10~ 7 


-9.74218 x 10" 8 
-9.74218 x 10" 8 




ypi [GeV] 


ReV ; K [GcV~ 4 ] 


lmV 0;K [GeV" 4 ] 


Our 
DK 


100 


-2.85709 x 10~ 7 
-2.85709 x 10~ 7 






Our 
DK 


200 


-7.61695 x 10" 8 
-7.61695 x 10~ 8 






Our 
DK 


300 


-3.29938 x 10~ 8 
-3.29938 x 10~ 8 






Our 
DK 


400 


-1.74228 x 10~ 8 
-1.74228 x 10~ 8 







Table 14: Comparison with the results of |12j (DK) in the setup of Ea. 1)3180 . Only the infrared finite part is shown. 
Unless indicated our relative error is below 10~ 5 . 
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s [{jex \ 




1m Ko;k l^ev J 


Uur 


1 A 1 

lU.z 




— lo.U5y(4J 






iu.o4oy 


— lo.uoyu 


Our 


y.z 


iy.yzo(o ) 


OO 1 TKf A \ 

— 11. Y f 0(4 ) 






m m on 


oo i7r;c; 
— zz.l (00 


Our 


o.z 


OA QC\QfQ\ 

z4.oyo(o ) 


— ZY.yoU(4 ) 


Do 




o/i nm £ 


—Li .y < oo 


UUI 


< .z 


oZ.loo(o ) 


Oft ^A1(0\ 

— oD.04 t[o) 


DS 




32.1805 


-36.5550 


Our 


6.2 


43.51(1) 


-50.114(8) 


DS 




43.4927 


-50.1010 


Our 


5.2 


62.62(2) 


-73.51(2) 


DS 




62.6575 


-73.5359 


Our 


4.2 


99.58(3) 


-120.09(2) 


DS 




99.6039 


-120.086 


Our 


3.2 


188.04(5) 


-237.07(5) 


DS 




188.017 


-237.028 



Table 15: Comparison with the results of |13| (DK) in the setup of Eg. 13191 . Only the infrared finite part is shown. 



Tensor coeff. 


Real part 


Imaginary part 


r c;K 


6.8986(3) x 10~ 6 


-6.3053(4) x 10~ 6 


V c;k\11 


-4.210(2) x 10~ 6 


3.423(3) x 10~ 6 


V c;k\12 


-2.9453(4) x 10~ 6 


2.4067(6) x 10" 6 


V° 

k c;k|11|11 


3.043(2) x 10~ 6 


-2.372(3) x 10~ 6 


V° 


2.1699(3) x 10~ 6 


-1.7005(4) x 10~ 6 


v c;k\12\12 


1.7392(2) x 10~ 6 


-1.5014(2) x 10" 6 



Table 16: Tensor Diagrams. Results for the tensor coefficients of rank 0, 1 and 2 for the diagram V c -k, see Ea. (|32()|l . 
with setup rrii = 7712 = mz — mi,, M — M z , m = ra e and s — M z (see Fig. 1101 . The unit of mass /1 is 1 GeV. Only 
the infrared finite part is shown. The results are in GeV -4 . 
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